1 Introduction

The question of universality, that is the independence of the critical properties of macroscopic systems from the microscopic details of the underlying model Hamiltonian, is a central issue in statistical physics, whose mathematical understanding is largely incomplete. A convenient framework where it can be studied is that of planar dimer models, which exhibit a rich critical behavior: algebraic decay of correlations, conformal invariance, and so on. The dimer model on a bipartite planar lattice is integrable and, more precisely, determinantal (also called ‘free fermionic’): its correlation functions are given by suitable minors of the so-called inverse Kasteleyn matrix [32]. The model is parametrized by edge weights \({\underline{t}}\) and has a non-trivial phase diagram. By varying \({\underline{t}}\), one can impose an average non-zero tilt \(\rho \) for the height field. A central object of the dimer model is the so-called characteristic polynomial P(zw), where zw are complex variables. For instance, the infinite-volume free energy is given by an integral of \(\log |P(z,w)|\) over the torus \({\mathbb {T}}=\{|z|=|w|=1\}\). Also, the large-distance decay of correlations is dictated by the so-called spectral curve, i.e. the algebraic curve \({\mathcal {C}}(P)=\{(z,w)\in {\mathbb {C}}^2:P(z,w)=0\}\). When the edge weights are such that the spectral curve intersects \({\mathbb {T}}\) transversally one is in the “liquid” or “massless” phase, where the two-point dimer–dimer correlation of the model decays like the inverse distance squared. Correspondingly the height field scales to a Gaussian Free Field (GFF) and the variance grows like the logarithm of the distance times \(1/\pi ^2\). Remarkably, this pre-factor is independent of the weights \({\underline{t}}\) and of the specific choice of the bipartite periodic planar lattice. This is related [34] to the fact that \({\mathcal {C}}(P)\) is a so-called Harnack curve. Summarizing, in the massless phase the scaling limit of height fluctuations of the dimer model is universal, in a very strong sense: the limit is always Gaussian, with logarithmic growth of the variance; moreover, the pre-factor in front of the logarithm in the variance is independent of the details of the underlying microscopic structure (edge weights and lattice).

The previous results heavily rely on the determinantal structure of the model, but universality is believed to hold much more generally. Motivated by this, we consider weak, translation-invariant, perturbations of the dimer model (for simplicity, we restrict to the square lattice). Generically, as soon as we switch on the perturbation, the determinantal structure provided by Kasteleyn’s theory breaks down. Two particular examples of perturbed, non-determinantal, models that we consider are: the 6-vertex model with general weights \(a_1,\dots , a_6\), in the disordered phase, close to, but not exactly at, the free-fermion point; and the dimer model with plaquette interaction, originally introduced in [29] and recently reconsidered in [2, 3, 38] in the context of quantum dimer models. There is a basic difference between these two cases: the 6-vertex model, even if non-determinantal, is still solvable via Bethe Ansatz (BA), see [4] and reference therein (note that the BA solution is not as explicit as the Kasteleyn solution of standard dimers: only a few thermodynamic functions can be explicitly computed). On the other hand, dimers with plaquette interaction are believed not to be solvable, i.e., not even the basic thermodynamic functions admit an explicit representation. From the exact solution, one finds that some of the critical exponents of the 6-vertex model depend continuously on the vertex weightsFootnote 1; they differ, in general, from those of the standard dimer model. On the other hand, the existence of non-trivial critical exponents in the dimer model with plaquette interaction, as well as in other planar models in the same ‘universality class’ (such as coupled Ising models, Ashkin–Teller and 8-vertex models) can be proved by constructive Renormalization Group (RG) methods [6, 9, 23, 24, 35], which allow one to express them as convergent power series in the interaction strength.

In this setting, it is natural to ask whether the height fluctuations are still described by a GFF at large scales and, in case, whether the pre-factor in front of the logarithm still displays some universal features. The very fact that the critical exponents depend non-trivially on the interaction strength suggests that universality cannot then be true in the naive, strong, sense that ‘large-scale properties are independent of the microscopic details of the model’: in fact, in this context, a weaker form of universality is expected, in the form of a number of scaling relations, originally proposed by Kadanoff [31], which allow one to determine all the critical exponents of the critical theory in terms of just one of them; this form of universality is often referred to as ‘weak universality’, see e.g. [4, Section 10.12]. Support for the Kadanoff scaling relations comes from the so-called bosonization picture, see e.g. [26] for a basic introduction. Only some of these universality relation have been rigorously proven [6, 9]; an example is the identity \(X_c X_e=1\) [31, Eq.(13b)], relating the “crossover exponent” \(X_c\) and “energy exponent” \(X_e\), see [6, Eq.(1.10)]. The proof in [6] covers both solvable and non-solvable models, but only works for scaling relations involving the critical exponents of the “local observables”, i.e., those that admit a representation in terms of a local fermionic operator. Other scaling relations, involving the critical exponents of non-local observables (e.g. monomer-monomer correlations in dimer models, or spin-spin correlations in the Ashkin–Teller model) remained elusive for many years. In particular, the relation \(X_p=X_e/4\) [31, Eq.(13a)], relating the energy exponent \(X_e\) to the “polarization exponent” \(X_p\) in the AT model, remains unproven at a rigorous level.

In this paper, we prove the stability of the Gaussian nature of the height fluctuations for non-integrable perturbations of the dimer model, with logarithmic growth of the variance in the whole liquid region. The pre-factor A in front of the logarithm depends, in general, non-trivially on the strength of the perturbation (see Remark 4 below) and on the dimer weights, so it is not universal in a naive, strong, sense. The non-trivial dependence of A on the interface tilt has been also verified numerically for the 6-vertex model [28]. Nevertheless, A satisfies a scaling relation, that connects it with the critical exponent of the dimer–dimer correlations.

Main Theorem

In a weakly perturbed dimer model with perturbation of strength \(\lambda \), the variance of the height difference between two faraway points grows like the logarithm of the distance, with a pre-factor \(A/\pi ^2\), where \(A=1+O(\lambda )\) is an analytic function of \(\lambda \) and of the dimer weights. Moreover, the prefactor satisfies the scaling relation

$$\begin{aligned} A=\nu , \end{aligned}$$
(1.1)

where \(2\nu \) is the anomalous decay exponent of the dimer–dimer correlation. Higher cumulants of the height difference between two points are bounded uniformly in their distance, that is, the fluctuations of the height difference are asymptotically Gaussian.

For a more precise statement, see Theorem 2 and the remarks and comments that follow it. Note that in the unperturbed case, \(\lambda =0\), the dimer–dimer correlation decays at large distances like \((dist.)^{-2}\) in the whole liquid phase, i.e., its decay exponent is equal to 2 (so that \(\nu =1\)), irrespective of the specific choice of the dimer weights. In this case, of course, our result reduces to the one of [34], \(A=1\). Note also that our result covers both integrable models, such as 6-vertex, and non-integrable ones, in the spirit of the universality picture.

Scaling relations involving exponents and amplitudes were conjectured by Haldane [30] and proved by Benfatto and Mastropietro [12, 13] in the context of quantum one-dimensional models. Even if formulated in different notations, the scaling relation (1.1) is strictly related to one of those proposed by Kadanoff, in particular to the above-mentioned, elusive, identity \(X_p=X_e/4\) [31, Eq.(13a)]. In fact, there is a duality (called ‘discrete bosonization’ in [18]) between the 6-vertex model, which is part of the class of perturbed dimer models considered in this paper, and the AT model; the duality implies non-trivial identities between the correlations of 6-vertex model and those of AT, see [18, Section 2.6]. In particular, the two-point correlation of the polarization operator in AT equals the ‘electric correlator’ \({\langle e^{i\pi (h_x-h_y)}\rangle }_{6V}\) of the 6-vertex model, see [18, Section 2.6]Footnote 2, while the energy critical exponent of AT equals the anomalous decay exponent of the arrow–arrow correlations of 6-vertexFootnote 3. Given these identities, (1.1) implies that \(X_p=X_E/4\) [31, Eq.(13a)], provided that

$$\begin{aligned} {\langle e^{i\pi (h_x-h_y)}\rangle }_{6V}\sim e^{-\frac{\pi ^2}{2}{\langle (h_x-h_y)^2\rangle }_{6V}}\sim e^{-\frac{A}{2}\log |x-y|} \end{aligned}$$
(1.2)

at large distances, as suggested by the asymptotic Gaussian behavior of the height differenceFootnote 4.

To prove our results, we start by periodizing the non-integrable dimer model on the toroidal graph of size L. Then we map it into a system of interacting two-dimensional lattice fermions, by rewriting its moment generating function as an integral over Grassmann variables, with non-quadratic action. At this point, we apply tools from the so-called constructive fermionic RG to control the \(L\rightarrow \infty \) limit of the correlation functions. In particular, we need a very sharp asymptotic description of the large-distance behavior of the dimer–dimer correlation function (cf. Theorem 1). The large-scale logarithmic behavior of height correlations, as well as the validity of the ‘Haldane’ scaling relation (1.1), rely on non-trivial identities (cf. (2.46)) between the coefficients appearing in the large-distance asymptotics of the dimer–dimer correlation function. In turn, (2.46) is the result of so-called Ward identities, i.e. exact relations between the correlation functions of the interacting lattice fermionic model, which the dimer model maps into.

The analogs of Theorems 1 and 2 have been proven in our previous works [25, 26] for the specific case of plaquette interaction and uniform edge weights \({\underline{t}}\equiv 1\). In this case, the average tilt of the height field is just \(\rho =0\) and the model has all the discrete symmetries of the lattice \({\mathbb {Z}}^2\). The extension to the general case, achieved here, is non-trivial: the loss of discrete rotation and reflection symmetries results, in the RG language, in the emergence of four new running coupling constants (two “Fermi velocities” and two “Fermi points”), whose flow, along the multi-scale integration procedure, has to be controlled via the choice of suitable counter-terms. Another consequence of the loss of rotation and reflection symmetry is that the cancellation at the basis of the logarithmic growth of the variance does not follow simply from the basic symmetries of the model, as it was the case in [25, 26]: the proof of the key identity, (2.46), now requires the use of a lattice Ward Identity for the dimer model, in combination with an emergent Ward Identity for an effective continuum model, which plays the role of ‘infrared fixed point’ of the RG flow. Quite surprisingly, the loss of rotation and reflection symmetry plays a role also in the technical control of the thermodynamic limit of correlations: in [25, 26], in order to simplify the analysis of the finite-size corrections to the critical correlation functions, we first studied a modified, slightly massive model of mass \(m>0\) (the modification consisted in adding a modulation of size m on the horizontal dimer weights; in the tilt-less case, this was enough to guarantee that the modified correlations decayed exponentially with rate m), and then we took the massless limit \(m\rightarrow 0\) after the thermodynamic limit. However, this strategy fails for general dimer weights: in this case, neither a modulation of the dimer weights nor other simple modifications of the model produce a mass; therefore, in the present paper, we directly derive quantitative estimates on the corrections to the thermodynamic limit of the massless correlations, by a careful control of the finite-size effects in the multi-scale procedure.

1.1 Related works

Let us conclude this introduction by mentioning some recent related works. While most literature on dimer models focuses on the determinantal case, there have been recently various suggestions to go beyond the exactly solvable situation [37]. As far as “limit shape phenomena” (i.e. laws of large numbers for the height profile) for non-solvable random interface models are concerned, let us mention for instance [15, 17, 36]. Closer in spirit to our results is [16], which provides a central limit theorem for height fluctuations of \(\nabla \phi \)-interface models with continuous heights and strictly convex potential. This work uses the Helffer–Sjöstrand formula, that is not available for discrete-height model like the dimer model. Let us mention also [1], which obtains convergence to the GFF (with interaction-dependent amplitude) for a dimer model with a special non-local interaction that makes it integrable, although not determinantal. Finally, a very interesting recent development is [14]: while in this work the convergence to the GFF is proven only for the non-interacting dimer model, the method of proof, that goes through Temperley’s bijection and Wilson’s algorithm rather than via Kasteleyn’s theory, might prove robust enough to allow for extensions to some non-determinantal situations.

1.2 Organization of the article

The rest of this work is organized as follows. The dimer model is defined in Sect. 2. There, we recall the large-scale behavior of the integrable model and we state our results for the non-integrable one. In Sect. 3 we give the Grassmann representation of the interacting dimer model and its lattice Ward identities. In Sect. 4 we recall the continuum reference model that plays the role of infrared fixed point of interacting dimers. Theorems 12 are proven in Sect. 5, conditionally on technical results, based on the multi-scale expansion, whose proofs are postponed to Sect. 6.

2 Model and Main Results

2.1 Dimers and height function

A dimer covering, or perfect matching, of a graph \(\Gamma \) is a subset of edges that covers every vertex exactly once. The set of dimer coverings of \(\Gamma \) is denoted \(\Omega _\Gamma \). We color the vertices of the bipartite graph \({\mathbb {Z}}^2\) black and white so that neighboring vertices have different colors. A white vertex is assigned the same coordinates \(x=(x_1,x_2)\) as the black vertex just at its left. The choice of coordinates is such that the vector \({{e}}_1\) is the one of length \(\sqrt{2}\) and angle \(-\pi /4\) w.r.t the horizontal axis, while \({{e}}_2\) is the one of length \(\sqrt{2}\) and angle \(+\pi /4\). The finite graph \({\mathbb {T}}_L\) denotes \({\mathbb {Z}}^2\) periodized (with period L) in both directions \({{e}}_1,{{e}}_2\). See Fig. 1.

Fig. 1
figure 1

The graph \({\mathbb {T}}_L\) for \(L=4\). The coordinate axes \({{e}}_1,{{e}}_2\), as well as the corresponding coordinates of some black/white vertices, are explicitly indicated.In the right drawing, the weights \(t_1,\dots ,t_4\) and thecorresponding edges of types \(1,\dots ,4\)

For simplicity we assume that L is even. Black/white sites are therefore indexed by coordinates \(x\in \Lambda =\{(x_1,x_2), 1\le x_i\le L\}\). An edge \(e=(b,w)\) of \({\mathbb {T}}_L\) is said to be of type \(r\in \{1,2,3,4\}\) if its white endpoint w is to the right, above, to the left or below the black endpoint b. If \(e=(b,w)\) is an edge of type r and x(b) is the coordinate of b then \(x(w)=x+v_r\), with

$$\begin{aligned} v_1=(0,0)\quad v_2=(-1,0)\quad v_3=(-1,-1)\quad v_4=(0,-1). \end{aligned}$$
(2.1)

If \(\Gamma \) is planar and bipartite, the height function allows us to interpret a dimer covering as a two-dimensional discrete surface. Let us recall the standard definition of height function for the infinite lattice \({\mathbb {Z}}^2\). Given \(M\in \Omega _{{\mathbb {Z}}^2}\), the height function \(h(\cdot ):=h_M(\cdot ) \) is defined on the dual lattice \(({\mathbb {Z}}^*)^2\), i.e. on the faces \(\eta \) of \({\mathbb {Z}}^2\). We set \(h(\eta _0):=0\) at a given reference face \(\eta _0\), and we let its gradients be given by

$$\begin{aligned} h(\eta ')-h(\eta )=\sum _{e\in C_{\eta \rightarrow \eta '}}\sigma _e ({\mathbb {1}}_{e}-1/4) \end{aligned}$$
(2.2)

where \(\eta ,\eta '\) are any two faces, \({\mathbb {1}}_{e}\) denotes the dimer occupancy, i.e., the indicator function that e is occupied by a dimer in M, while \(C_{\eta \rightarrow \eta '}\) is any nearest-neighbor path on the dual lattice \(({\mathbb {Z}}^*)^2\) from \(\eta \) to \(\eta '\) (the right side of (2.2) is independent of the choice of \(C_{\eta \rightarrow \eta '}\)). The sum runs over the edges crossed by the path and \(\sigma _e=+1/-1\) depending on whether the oriented path \(C_{\eta \rightarrow \eta '}\) crosses e with the white site on the right/left.

2.2 Definition of the model

We define here both the non-interacting dimer model [32] and the interacting one. Both are probability measures on \(\Omega _L:=\Omega _{{\mathbb {T}}_L}\), denoted \({\mathbb {P}}_{L,{\underline{t}}}\) and \({\mathbb {P}}_{L,\lambda , {\underline{t}}}\) respectively, where \(\lambda \in {\mathbb {R}}\) is the interaction strength and \({\underline{t}}\) are the edge weights. For lightness of notation, the index \({\underline{t}}\) will be dropped.

2.2.1 The non-interacting dimer model

We assign a positive weight to each edge. More precisely, an edge of type \(r\in \{1,2,3,4\}\) is given a weight \(t_r> 0\). Then, the weight of a configuration \(M\in \Omega _L\) is

$$\begin{aligned} {\mathbb {P}}_L(M)= & {} \frac{ t_1^{N_1(M)}t_2^{N_2(M)}t_3^{N_3(M)}t_4^{N_4(M)}}{Z^0_L}, \end{aligned}$$
(2.3)
$$\begin{aligned} Z^0_L= & {} \sum _{M'\in \Omega _L} t_1^{N_1(M')}t_2^{N_2(M')}t_3^{N_3(M')}t_4^{N_4(M')} \end{aligned}$$
(2.4)

with \(N_i(M)\) the number of dimers on edges of type i in configuration M. Since the total number of dimers is constant, we can rescale all weights by a common factor and we will set \(t_4\equiv 1\) from now on. It is known that the free energy per site has a limit as \(L\rightarrow \infty \) (the infinite volume free energy):

$$\begin{aligned} F({\underline{t}})= & {} \lim _{L\rightarrow \infty }\frac{1}{L^2} \log Z^0_L=\frac{1}{(2\pi )^2}\int _{[-\pi ,\pi ]^2}d k \log \mu (k), \end{aligned}$$
(2.5)
$$\begin{aligned} \mu (k)= & {} t_1+ it_2 e^{i k_1}-t_3e^{i k_1+i k_2}-i e^{i k_2}. \end{aligned}$$
(2.6)

Note that

$$\begin{aligned} \mu (k)=\mu ^*((\pi ,\pi )-k). \end{aligned}$$
(2.7)

The “characteristic polynomial” mentioned in the introduction is \(P(z,w):=\mu (-i \log z, -i \log w)\).

Also, the measure \({\mathbb {P}}_L\) itself has a limit \({\mathbb {P}}\) as \(L\rightarrow \infty \), in the sense that the probability of any local event converges. The non-interacting model is integrable, and both the measure \({\mathbb {P}}_L\) and its limit \({\mathbb {P}}\) admit a determinantal representation, recalled in Sect. 3.1.

In the special case where \(t_1=t_3=: t\) and \(t_2=1\), i.e. assigning weight t to horizontal edges and 1 to vertical ones, one recovers the model originally solved by Kasteleyn [32]. For general weights \(t_1,t_2,t_3\), the model is equivalent to Kasteleyn’s model with different weights for horizontal and vertical edges, and a non-zero average slope \(\rho =\rho (t_1,t_2,t_3)\in {\mathbb {R}}^2\) for the height function, i.e.,

$$\begin{aligned} {\mathbb {E}}(h(\eta +{{e}}_i)-h(\eta ))=\rho _i,\qquad i=1,2, \end{aligned}$$
(2.8)

where \({\mathbb {E}}\) denotes the average with respect to \({\mathbb {P}}\). In fact, the weights \(t_i\) are chemical potentials by which one can fix the densities of the four types of edges. Then, the slope \(\rho \) is obtained as a function of the four densities using the definition (2.2) of height function.

Another special case is obtained letting e.g. \(t_3\rightarrow 0\): then, the model reduces to the closed-packed dimer model on the hexagonal graph with weights \(1,t_1,t_2\) for the three types of edges.

Note that the condition \(\mu (k)=0\) gives

$$\begin{aligned} e^{i k_2}=\frac{t_1+i t_2 e^{i k_1}}{i+ t_3 e^{i k_1}} \end{aligned}$$
(2.9)

that determines the intersections of two circles in the complex plane. We will make the following important assumption:

Assumption 1

The parameters \({\underline{t}}\) are such that \(\mu (\cdot )\) has two distinct simple zeros, that we call \(p^+\) and \(p^-\), on \([-\pi ,\pi ]^2\) (i.e. the two circles intersect transversally). In view of (2.7), one has \(p^++p^-=(\pi ,\pi )\).

Remark 1

Note that, under Assumption 1, none of the weights \(t_1,t_2,t_3,1\) exceeds the sum of the other three, otherwise \(\mu (k)\) would vanish nowhere on \([-\pi ,\pi ]^2\). Note also that \(p^\omega ,\omega =\pm \) cannot coincide with any of the four values \(k=(\epsilon _1\pi /2,\epsilon _2\pi /2), \epsilon _1=\pm 1, \epsilon _2=\pm 1\), otherwise one would have \(p^+=p^-\) (modulo \((2\pi ,2\pi )\)).

Under Assumption 1, it is known [34] that the infinite-volume measure has power-law decaying correlations (in the language of [34], the dimer model is said to be in a “liquid phase”). With the nomenclature of condensed matter theory, the zeros \(p^\pm \) are called “Fermi points”.

2.3 The interacting dimer model, and relation to the 6-vertex model

In order to study the effect of the breaking of integrability we introduce interacting dimer measures of the following form:

$$\begin{aligned} {\mathbb {P}}_{L,\lambda }(M)= \frac{p_{L,\lambda }(M)}{Z_{L}} \end{aligned}$$
(2.10)

where

$$\begin{aligned} \begin{aligned}&p_{L,\lambda }(M)=t_1^{N_1(M)}t_2^{N_2(M)}t_3^{N_3(M)} \, e^{\lambda W_L(M)},\\&Z_{L}=\sum _{M\in \Omega _L}p_{L,\lambda }(M) \end{aligned} \end{aligned}$$
(2.11)

and the interaction potential \(W_L\) is given as

$$\begin{aligned} W_L(M)=\sum _{x\in \Lambda }f(\tau _x M), \end{aligned}$$
(2.12)

where f is some fixed local function of the dimer configuration and \(\tau _x M\) denotes the configuration M translated by \(x_1 {{e}}_1+x_2{{e}}_2\). We do not require \(f(\cdot )\) to be symmetric under reflections or rotation by \(\pi /4\).

Let us mention two interesting particular examples of interaction \(W_L(M)\). The first one is the plaquette interaction that was considered in our works [25, 26] and previously in the theoretical physics literature [2] in the context of quantum dimer models. Namely,

$$\begin{aligned} W_L(M)=\sum _{\eta \in {\mathbb {T}}^*_L}{} \mathbf{1}_\eta (M) \end{aligned}$$
(2.13)

where the sum runs over all faces of \({\mathbb {T}}_L\) and \(\mathbf{1}_\eta (M)\) is the indicator function that two of the four edges surrounding \(\eta \) are occupied by dimers. In this case the function f in (2.12) is

$$\begin{aligned} f_P(M)={\mathbb {1}}_{e_1}{\mathbb {1}}_{e_2} +{\mathbb {1}}_{e_3}{\mathbb {1}}_{e_4} +{\mathbb {1}}_{e_1}{\mathbb {1}}_{e_5}+{\mathbb {1}}_{e_6}{\mathbb {1}}_{e_7} \end{aligned}$$
(2.14)

with \(e_1,\dots ,e_7\) as in Fig. 2.

Fig. 2
figure 2

The edges appearing in (2.14). b denotes the black vertex of coordinates (0, 0)

Another important example (see again Fig. 2) is

$$\begin{aligned} f_{6v}(M):={\mathbb {1}}_{e_1}{\mathbb {1}}_{e_2} +{\mathbb {1}}_{e_3}{\mathbb {1}}_{e_4}. \end{aligned}$$
(2.15)

In this case, the interaction \(W_L(M)\) in (2.13) is modified in that the sum runs only over one of the two sub-lattices of \({\mathbb {T}}_L^*\) (the subset of faces with black top-right vertex). Then, it is known that this interacting dimer model is equivalent to the 6-vertex model [5, 19, 20]. Recall that configurations of the 6-vertex model are assignments of orientations (arrows) to the edges of \({\mathbb {Z}}^2\) such that at each vertex there are two incoming and two outgoing arrows. There are 6 possible arrow configurations at any vertex, each being assigned a positive weight \(a_1,\dots ,a_6\) (see Fig. 3) and the weight of a configuration is the product of the weights over all vertices.

Fig. 3
figure 3

The six possible vertex configurations of the 6-vertex model and the associated weights

By multiplying all weights by a common factor, one can reduce e.g. to \(a_3=1\). Moreover, on the torus, the number of vertices of type 5 equals the number of vertices of type 6, so one can set without loss of generality \(a_5=1\). One is left with four positive weights \(a_1,a_2,a_4,a_6\) and the model can be mapped to the interacting dimer model with weights \(t_1,t_2,t_3\), interaction (2.15) and interaction parameter \(\lambda \) such that

$$\begin{aligned} t_1=a_1,\;t_2=a_4,\;t_3=a_2, \;(t_1 t_3+t_2)e^\lambda =a_6. \end{aligned}$$
(2.16)

More precisely, as in Fig. 4, the dimer model lives on a square grid rotated by 45 degrees w.r.t. the lattice of the 6-vertex model.

Fig. 4
figure 4

The 6-vertex model lives on the square grid \({\mathcal {G}}_{6v}\) with dotted edges, while the dimer model lives on the square grid \({\mathcal {G}}_{d}\) with full edges. Faces of \({\mathcal {G}}_{d}\) containing a vertex of \({\mathcal {G}}_{6v}\) are called “even faces” and the others “odd faces”

The mapping is obtained by associating to the arrow configuration at a vertex x of \({\mathcal {G}}_{6v}\) a dimer configuration at the even face of \({\mathcal {G}}_{d}\) containing x, as in Fig. 5.

Fig. 5
figure 5

The local arrow-to-dimer mapping

The map is one-to-many because arrow configurations of type 6 are mapped to two possible dimer configurations. However, it is easily checked that the partition functions of the two models are equal provided the parameters are identified as in (2.16). Moreover, the height function of the dimer model, restricted to odd faces of \({\mathcal {G}}_d\), equals (up to a global prefactor) the canonical height function of the 6-vertex model [40]. The 6-vertex model is known to be free-fermionic (i.e. determinantal) if and only if

$$\begin{aligned} \Delta :=\frac{a_1 a_2+a_3 a_4-a_5a_6}{2\sqrt{a_1 a_2 a_3 a_4}}=0. \end{aligned}$$
(2.17)

It is immediately checked that this condition is equivalent to \(\lambda =0\) for the interacting dimer model.

2.4 Non-interacting model: dimer–dimer correlations and logarithmic height fluctuations

It is known [34] that, under the infinite-volume measure \({\mathbb {P}}\) of the non-interacting model, dimer–dimer correlations decay like the inverse distance squared and the height field behaves on large scales like a massless Gaussian field. We briefly recall the basic facts here, since they serve to motivate our main result for the interacting dimer model. For \(\omega =\pm \), we let

$$\begin{aligned} \alpha _\omega&=\partial _{k_1}\mu (p^\omega )=-t_2e^{i p^\omega _1}-i t_3 e^{i(p^\omega _1+p^\omega _2)}=-i t_1- e^{i p^\omega _2}, \end{aligned}$$
(2.18)
$$\begin{aligned} \beta _\omega&=\partial _{k_2}\mu (p^\omega )=-i t_3 e^{i (p^\omega _1+p^\omega _2)}+e^{i p^\omega _2}=-i t_1+t_2e^{i p^\omega _1}, \end{aligned}$$
(2.19)

where \(p^\pm \) are the two zeros of \(\mu (\cdot )\), as in Assumption 1. (The complex numbers \(\alpha _\omega ,\beta _\omega \) are called “Fermi velocities” in the jargon of condensed matter.) Define also

$$\begin{aligned} \phi _\omega :x\in {\mathbb {R}}^2\mapsto \phi _\omega (x):=\omega (\beta _\omega x_1 -\alpha _\omega x_2)\in {\mathbb {C}}. \end{aligned}$$
(2.20)

Remark 2

Under Assumption 1 on the weights \({\underline{t}}\), the complex numbers \(\alpha _\omega \) and \(\beta _\omega \) are not colinear, as elements of the complex plane [34], i.e. \(\alpha _\omega /\beta _\omega \) is not real. Therefore, \(\phi _\omega \) is a bijection from \({\mathbb {R}}^2\) to the complex plane. More precisely, one has that

$$\begin{aligned} {{\mathrm{Im}}}(\beta _+/\alpha _+)>0. \end{aligned}$$
(2.21)

In fact, parametrize the weights \(t_1,t_2,t_3\) as

$$\begin{aligned} t_1=t e^{-B_1},\quad t_2=e^{-B_1-B_2},\quad t_3=t e^{-B_2},\quad B_1,B_2\in {\mathbb {R}}. \end{aligned}$$

For \(B_1,B_2=0\) it is immediately checked that \(p^+=(0,0),p^-=(\pi ,\pi )\) and that (2.21) holds. On the other hand, once t is fixed, it is known [34] that the set \({\mathcal {B}}_t\) of values of \(B=(B_1,B_2)\) for which Assumption 1 holds is a connected subset of \({\mathbb {R}}^2\) on which \({\mathrm{Im}}(\beta _+/\alpha _+)\) vanishes nowhere, and it is therefore everywhere positive.

Because of the symmetry (2.7), one has \(\alpha _\omega =-\alpha ^*_{-\omega }\), \(\beta _{\omega }=-\beta ^{*}_{-\omega }\) and \(\phi _\omega ^*(\cdot )=\phi _{-\omega }(\cdot )\).

The relation between the massless Gaussian field and the height function is given by the following results. Let n be an integer and \(\eta _j,j\le 2n\) be faces of \({\mathbb {Z}}^2\). With some abuse of notation, we identify a face \(\eta \) with its mid-point. Then,

$$\begin{aligned}&{\mathbb {E}}\left[ (h(\eta _1)-h(\eta _2));(h(\eta _3)-h(\eta _4)) \right] \nonumber \\&\quad =\frac{1}{2\pi ^2} \mathfrak {R}\log \left( \frac{(\phi _+(\eta _4)-\phi _+(\eta _1))(\phi _+(\eta _3) - \phi _+(\eta _2))}{(\phi _+(\eta _4)-\phi _+(\eta _2))(\phi _+(\eta _3)-\phi _+(\eta _1))} \right) \nonumber \\&\qquad +O\left( \frac{1}{\min _{i\ne j\le 4}|\eta _i-\eta _j|+1}\right) \end{aligned}$$
(2.22)

where \(\phi _+(\eta _i)-\phi _+(\eta _j)\) should be read as 1 in case \(\eta _i=\eta _j\). Also, for \(n>2\)

$$\begin{aligned}&{\mathbb {E}}\left[ (h(\eta _1)-h(\eta _2));\dots ; (h(\eta _{2n-1})-h(\eta _{2n})) \right] \nonumber \\&\quad =O\left( \frac{1}{\min _{i\ne j\le 2n}|\eta _i-\eta _j|+1}\right) \end{aligned}$$
(2.23)

where \({\mathbb {E}}(X_1;\dots ;X_k)\) denotes the joint cumulant of the random variables \(X_1,\dots ,X_k\). In particular, as \({|\eta _1-\eta _2|\rightarrow \infty }\),

$$\begin{aligned} {\mathrm{Var}}_{{\mathbb {P}}}(h(\eta _1)-h(\eta _2))=\frac{1}{\pi ^2}\mathfrak {R}\log (\phi _+(\eta _1)-\phi _+(\eta _2))+O(1) \end{aligned}$$
(2.24)

while the cumulants of order \(n\ge 3\) of \((h(\eta _1)-h(\eta _2))\) are bounded from above, uniformly in \(\eta _1,\eta _2\). It is well known that (2.22) and (2.23) imply that the height field tends, in the scaling limit, to a GFF with covariance

$$\begin{aligned} -\frac{1}{2\pi ^2}\mathfrak {R}\log (\phi _+({ x})-\phi _+({y})), \end{aligned}$$
(2.25)

with \(x,y\in {\mathbb {R}}^2\). For (2.22) see [34] and for (2.23) see e.g. [25, Th. 5] (in [25] the weights \(t_i\) are all 1 and \(\eta _1=\eta _3=\dots \eta _{2n-1}, \eta _2=\eta _4=\dots =\eta _{2n}\); the proof of (2.23) in the general case works the same way).

Remark 3

Note that the prefactor \(1 /\pi ^2\) is independent of the weights \({\underline{t}}\). In [34], such universality is related to the fact that the spectral curve, i.e. the algebraic curve defined by the zeros on \({\mathbb {C}}^2\) of the polynomial \(P(z,w):=\mu (-i \log z,-i\log w)\), is a so-called Harnack curve.

It is useful to recall the key points of the proof of (2.22) in order to understand the main new features posed by the presence of the interaction. From the definition of height function,

$$\begin{aligned}&{\mathbb {E}}\left[ (h(\eta _1)-h(\eta _2));(h(\eta _3)-h(\eta _4)) \right] \nonumber \\&\quad = \sum _{e\in C_{\eta _1\rightarrow \eta _2}} \sum _{e'\in C_{\eta _3\rightarrow \eta _4}}\sigma _e \sigma _{e'} {\mathbb {E}}({\mathbb {1}}_e;{\mathbb {1}}_{e'}) \end{aligned}$$
(2.26)

where \({\mathbb {E}}({\mathbb {1}}_e;{\mathbb {1}}_{e'})\) is the dimer–dimer correlation function

$$\begin{aligned} {\mathbb {E}}({\mathbb {1}}_e;{\mathbb {1}}_{e'}):= {\mathbb {E}}({\mathbb {1}}_e{\mathbb {1}}_{e'})- {\mathbb {E}}({\mathbb {1}}_e){\mathbb {E}}({\mathbb {1}}_{e'}). \end{aligned}$$

This correlation function has an exact expression involving the inverse Kasteleyn matrix of the infinite lattice; at large distances, it can be expressed as

$$\begin{aligned}&{\mathbb {E}}({\mathbb {1}}_e;{\mathbb {1}}_{e'})= A_{r,r'} (x,x')+B_{r,r'}(x,x')+R_{r,r'}(x,x'),\nonumber \\&A_{r,r'}(x,x')=\frac{1}{4\pi ^2}\sum _{\omega =\pm } \frac{K_{\omega , r} K_{\omega , r'}}{(\phi _\omega (x-x'))^2}\nonumber \\&B_{r,r'}(x,x')=\frac{1}{4\pi ^2} \sum _{\omega =\pm } \frac{K_{-\omega , r} K_{\omega ,r'} }{|\phi _\omega (x-x')|^2} e^{i (p^\omega -p^{-\omega })\cdot (x-x')} \end{aligned}$$
(2.27)

where:

  • the edge e (resp. \(e'\)) is of type \(r=r(e)\) (resp. \(r'=r(e')\)) and the coordinate of its black endpoint is \(x=x(e)\) (resp. \(x'=x(e')\));

  • \(K_{\omega ,r}= K_r e^{-i p^\omega \cdot v_r}\) (see (2.1) for the definition of \(v_i\)) with

    $$\begin{aligned} K_1=t_1,\quad K_2=i t_2, \quad K_3=-t_3, \quad K_4=-i; \end{aligned}$$
    (2.28)

    note that \(K_{-\omega ,r}=K^*_{\omega ,r}\).

  • \(R_{r,r'}(x,x')\) is a remainder, decaying like \(|x-x'|^{-3}\) at large distance.

By plugging (2.27) into (2.26), one obtains, letting \(d:=\min _{i\ne j}|\eta _i-\eta _j|\),

$$\begin{aligned} (2.26)= & {} \frac{1}{4\pi ^2} \sum _{\omega =\pm } \sum _{\begin{array}{c} e\in C_{\eta _1\rightarrow \eta _2}\\ e'\in C_{\eta _3\rightarrow \eta _4} \end{array}}\sigma _e \sigma _{e'} \Big [\frac{K_{\omega , r} K_{\omega , r'}}{(\phi _\omega (x-x'))^2}\nonumber \\&+\frac{K_{-\omega , r} K_{\omega ,r'} }{|\phi _\omega (x-x')|^2} e^{i (p^\omega -p^{-\omega })\cdot (x-x')}\Big ]+O\Big (\frac{1}{d+1}\Big ), \end{aligned}$$
(2.29)

where the \(O(\frac{1}{d+1})\) in the second line comes from the remainder \(R_{r,r'}(x,x')\) (in order to prove this error bound, one needs to choose the paths \(C_{\eta _1\rightarrow \eta _2}\), \(C_{\eta _3\rightarrow \eta _4}\) to be as ‘well separated’ as possible, cf. with [25, Fig.3] and the discussion in [25, Sect.3.2]).

In order to conveniently rewrite the double sum over \(e,e'\) in (2.29), we assume for simplicity that the paths \(C_{\eta _1\rightarrow \eta _2},C_{\eta _3\rightarrow \eta _4}\) are a concatenation of elementary steps in direction \(\pm {{e}}_1\) and \(\pm {{e}}_2\), connecting faces of the same parity: e.g., assume that an elementary step s(x, 1) in direction \(+{{e}}_1\) ‘centered at x’ consists in crossing the two bonds \((\bullet ,\circ )=(x,x+v_3)\) and \((\bullet ,\circ )=(x,x+v_4)\) with the white vertex on the right, while an elementary step s(x, 2) in direction \(+{{e}}_2\) centered at x consists in crossing the two bonds \((\bullet ,\circ )=(x,x)\) and \((\bullet ,\circ )=(x-v_4,x)\) with the white vertex on the right. A simple but crucial observation is that

$$\begin{aligned}&\sum _{e\in s(x,1)} \sigma _e K_{\omega ,r(e)} = K_{3}e^{-i p^\omega v_3} +K_{4}e^{-i p^\omega v_4} =-i \beta _\omega = -i \omega \Delta _1\phi _\omega \end{aligned}$$
(2.30)
$$\begin{aligned}&\quad \sum _{e\in s(x,2)} \sigma _e K_{\omega ,r(e)} = K_1e^{-ip^\omega v_1}+K_4e^{-ip^\omega v_4} =i\alpha _\omega = -i \omega \Delta _2\phi _\omega \end{aligned}$$
(2.31)

where \(\Delta _j\phi _\omega \) denotes the discrete gradient in direction \({{e}}_j\) of the linear function \(\phi _\omega \) defined in (2.20). Therefore, by using these identities, we can rewrite the right side of (2.29), up to the error \(O(\frac{1}{d+1})\), as

$$\begin{aligned}&-\frac{1}{4\pi ^2} \sum _{\omega =\pm } \sum _{\begin{array}{c} s(x,j)\in C_{\eta _1\rightarrow \eta _2}\\ s(x',j')\in C_{\eta _3\rightarrow \eta _4} \end{array}}\frac{\Delta _j\phi _\omega (x)\Delta _{j'}\phi _\omega (x')}{(\phi _\omega (x-x'))^2}\nonumber \\&+\frac{1}{4\pi ^2} \sum _{\omega =\pm } \sum _{\begin{array}{c} s(x,j)\in C_{\eta _1\rightarrow \eta _2}\\ s(x',j')\in C_{\eta _3\rightarrow \eta _4} \end{array}}\frac{\Delta _j\phi _{-\omega }(x)\Delta _{j'}\phi _\omega (x')}{|\phi _\omega (x-x')|^2} e^{i (p^\omega -p^{-\omega })\cdot (x-x')}. \end{aligned}$$
(2.32)

Note that, since \(p^\pm \) are distinct by assumption, the summand in the second line is genuinely oscillating. By using the oscillations, one finds that the contribution from the second line can be bounded by \(O(\frac{1}{d+1})\), like the one from the sum over \(R_{r,r'}(x,x')\); cf. with the discussion in [25, Sect.3.2]. We are left with the term in the first line, which is the Riemann sum approximation of the integral in the complex plane

$$\begin{aligned} - \frac{1}{2\pi ^2}\mathfrak {R}\int _{\phi _+(\eta _1)}^{\phi _+(\eta _2)}dz\int _{\phi _+(\eta _3)}^{\phi _+(\eta _4)}dz' \frac{1 }{(z-z')^2} \end{aligned}$$
(2.33)

whose explicit evaluation gives the main term in the r.h.s. of (2.22).

2.5 The interacting case: main results

In the presence of the interaction, \(\lambda \not =0\), Kasteleyn theory is not valid anymore, so that one cannot rely on an explicit computation of the dimer correlations to check the validity of the asymptotic Gaussian behavior of the height function. However, dimer correlations can be written as a renormalized expansion based on multiscale analysis. From now on, we will assume that the interaction is small:

$$\begin{aligned} |\lambda |\le \varepsilon \end{aligned}$$
(2.34)

and all claims above hold if \(\varepsilon \) is small enough (uniformly in L).

Our first result is:

Theorem 1

Given a local function g of the dimer configuration, the limit

$$\begin{aligned} {\mathbb {E}}_\lambda (g):=\lim _{L\rightarrow \infty } {\mathbb {E}}_{L,\lambda }(g) \end{aligned}$$
(2.35)

exists. The infinite-volume dimer–dimer correlations are given by

$$\begin{aligned}&{\mathbb {E}}_{\lambda }({\mathbb {1}}_{e};{\mathbb {1}}_{e'})= {\bar{A}}_{r,r'}(x,x')+{\bar{B}}_{r,r'}(x,x')+{\bar{R}}_{r,r'}(x,x') \end{aligned}$$
(2.36)
$$\begin{aligned}&{\bar{A}}_{r,r'}(x,x')=\frac{1}{4\pi ^2} \sum _{\omega =\pm } \frac{{\bar{K}}_{\omega , r} {\bar{K}}_{\omega , r'} }{{\bar{\phi }}_\omega (x-x')^2} \end{aligned}$$
(2.37)
$$\begin{aligned}&{\bar{B}}_{r,r'}(x,x')=\frac{1}{4\pi ^2} \sum _\omega \frac{{\bar{H}}_{-\omega , r} {\bar{H}}_{\omega , r'} }{|{\bar{\phi }}_\omega (x-x')|^{2\nu } } e^{i ({\bar{p}}^\omega -{\bar{p}}^{-\omega })\cdot (x-x')} \end{aligned}$$
(2.38)

where:

  • \(r=r(e)\) is the type of the edge e, \(x=x(e)\) is the coordinate of the black site of e, and similarly for \(r',x'\);

  • \({\bar{\phi }}_\omega (x)=\omega ({\bar{\beta }}_\omega x_1-{\bar{\alpha }}_\omega x_2)\);

  • one has

    $$\begin{aligned}&\nu =1+O(\lambda )\in {\mathbb {R}}, \end{aligned}$$
    (2.39)
    $$\begin{aligned}&{\bar{K}}_{\omega , r}=K_{\omega , r}+O(\lambda )\in {\mathbb {C}},\quad {\bar{H}}_{\omega , r}= K_{\omega , r}+O(\lambda )\in {\mathbb {C}} \end{aligned}$$
    (2.40)
    $$\begin{aligned}&{\bar{\alpha }}_\omega =\alpha _\omega +O(\lambda )\in {\mathbb {C}},\quad {\bar{\beta }}_\omega =\beta _\omega +O(\lambda )\in {\mathbb {C}}, \end{aligned}$$
    (2.41)
    $$\begin{aligned}&{\bar{p}}^\omega = p^\omega +O(\lambda )\in [-\pi ,\pi ]^2; \end{aligned}$$
    (2.42)

    these are all analytic functions of \(\lambda \) and satisfy the symmetries

    $$\begin{aligned}&{\bar{\alpha }}_\omega ^*=-{\bar{\alpha }}_{-\omega },\quad {\bar{\beta }}_\omega ^*=-{\bar{\beta }}_{-\omega }, \end{aligned}$$
    (2.43)
    $$\begin{aligned}&{\bar{K}}_{\omega , r}^*={\bar{K}}_{-\omega ,r},\quad {\bar{H}}_{\omega , r}^*={\bar{H}}_{-\omega ,r} \end{aligned}$$
    (2.44)
    $$\begin{aligned}&{\bar{p}}^++{\bar{p}}^-=(\pi ,\pi ). \end{aligned}$$
    (2.45)

Finally, \({\bar{R}}_{r,r'}(x,x')=O(|x-x'|^{-5/2})\) (the exponent 5/2 could be replaced by any \(\delta <3\) provided \(\lambda \) is small enough).

[A warning on notation: given a quantity (such as \(\alpha _\omega ,\phi _\omega \)) referring to the non-interacting model, the corresponding \(\lambda \)-dependent quantity for the interacting model will be distinguished by a bar, such as \({\bar{\alpha }}_\omega ,\) etc. On the other hand, we denote by \(z^*\) the complex conjugate of a number z.]

Note that the interaction modifies the decay rate of the correlation, producing a non-trivial (‘anomalous’) critical exponent \(\nu \). The analytic functions appearing in (2.39) are expressed as convergent power series but, due to the complexity of the expansion, the coefficients can be explicitly evaluated only at the lowest orders. This makes impossible to verify directly the validity of relations like (2.30), which were essential for the proof of large-scale Gaussian behavior of the height field in the non-interacting case. However, we can prove non-perturbatively that the parameters appearing in (2.36) are not independent, but related by exact relations, which are the central result of the present work:

Theorem 2

One has

$$\begin{aligned} \sum _{e\in s(x,j)} \sigma _e {\bar{K}}_{\omega ,r(e)}=-i\omega \sqrt{\nu }\Delta _j{\bar{\phi }}_\omega , \end{aligned}$$
(2.46)

where \(\nu =\nu (\lambda )\) is the same as the critical exponent in Theorem 1. Here, s(xj) is the elementary step in direction \(+{{e}}_{j}\) centered at x, thought of as a collection of two bonds, as defined before (2.30). As a consequence,

$$\begin{aligned}&{\mathbb {E}}_\lambda \left[ (h(\eta _1)-h(\eta _2));(h(\eta _3)-h(\eta _4)) \right] \nonumber \\&\quad =\frac{\nu }{2\pi ^2} \mathfrak {R}\log \left( \frac{({\bar{\phi }}_+(\eta _4)-{\bar{\phi }}_ +(\eta _1))({\bar{\phi }}_+(\eta _3) - {\bar{\phi }}_+(\eta _2))}{({\bar{\phi }}_+(\eta _4) -{\bar{\phi }}_+(\eta _2))({\bar{\phi }}_+(\eta _3)-{\bar{\phi }}_+(\eta _1))} \right) \nonumber \\&\qquad +O\left( \frac{1}{\min _{i\ne j\le 4}|\eta _i-\eta _j|^{1/2}+1}\right) \end{aligned}$$
(2.47)

(the exponent 1/2 could be replaced by any \(\delta <1\) provided \(\lambda \) is small enough; as in (2.22), when \(\eta _i=\eta _j\), \({\bar{\phi }}_+(\eta _i)-{\bar{\phi }}_+(\eta _j)\) has to be read as 1).

Note that the result contains two non-trivial pieces of information: first, the sum of \(\sigma _e {\bar{K}}_{\omega ,r(e)}\) along a step in direction \({{e}}_i\) is proportional to the discrete gradient of \({\bar{\phi }}_\omega \) in the same direction; second, the coefficient of proportionality is related in an elementary way to the critical exponent \(\nu \) that appears in (2.38). The latter relation immediately implies (see the argument spelled out after (2.26) for the non-interacting model, and in Section 5.1 below for the interacting case) the identity (cf. (2.47)) between height fluctuation amplitude and critical exponent \(\nu \) and is a form of universality.

Remark 4

Recall that for the non-interacting model \(\nu =1\), in particular it is independent of the weights \(t_i\). This is not true any more for the interacting model. Indeed, an explicit calculation of \(\nu \) at first order in \(\lambda \) for the model with plaquette interaction shows a non-trivial dependence both on \(\lambda \) and on the weights [27].

Theorem 2 follows from a combination of exact relations among correlation functions of the interacting dimer model (“lattice Ward identities”) together with chiral gauge symmetry emerging in the continuum scaling limit; it is remarkable that such a symmetry, valid only in the continuum limit, implies nevertheless exact relations for the coefficients of the lattice theory.

Remark 5

The analog of Theorem 1 has been proven in [25, 26] in the special case \(t_1=t_2=t_3=1\) and with plaquette interaction as in (2.14), which has the same discrete symmetries as the lattice. In that case, for symmetry reasons one obtains automatically that the ratios \( \frac{{\bar{K}}_{\omega ,r}}{K_{\omega ,r}} \) are independent of \(r,\omega \) and that \( \frac{{\bar{\alpha }}_\omega }{\alpha _\omega }=\frac{{\bar{\beta }}_\omega }{\beta _\omega }, \) the ratios being again \(\omega \)-independent.

Let us add also that, in the works [25, 26], the existence of the \(L\rightarrow \infty \) limit of the measure \({\mathbb {P}}_{L,\lambda }\) itself was not proven: instead, we modified the measure \({\mathbb {P}}_{L,\lambda }\) by an infra-red cut-off \(m>0\) (mass) and then we took the limit where first \(L\rightarrow \infty \) and then \(m\rightarrow 0\). We explain in Sect. 6 how the need of the cut-off m can be bypassed.

To upgrade Theorem 2 into a statement of convergence of the height field to a Gaussian Free Field with covariance

$$\begin{aligned} -\frac{\nu }{2\pi ^2}\mathfrak {R}\log ({\bar{\phi }}_+(x)-{\bar{\phi }}_+(y)), \end{aligned}$$

one needs to complement (2.47) with the statement that higher cumulants are negligible, i.e. that, for \(n>2\) and some \(\theta >0\),

$$\begin{aligned} {\mathbb {E}}_\lambda \left[ (h(\eta _1)-h(\eta _2));\ldots ;(h(\eta _{2n-1})-h(\eta _{2n})) \right] =O((\min _{i\ne j}|\eta _i-\eta _j|+1)^{-\theta }). \end{aligned}$$

In turn, this requires an analog of (2.36) for multi-dimer correlation functions. This can be done following the ideas of Sects. 5 and 6 below but, in order to keep this work within reasonable length, we decided not to develop this point. The interested reader may look at [25, Theorem 3 and Sec. 7], where the precise statements on multi-dimer correlations and on the convergence to the GFF are given in detail for the model with edge weights \({\underline{t}}\equiv 1\) and interaction (2.13).

3 Grassmann Integral Representation

3.1 Kasteleyn theory

For the statements of this section and more details on Kasteleyn theory, we refer the reader for instance to [33, 34].

The partition function and the correlations of the non-interacting model (2.3) can be explicitly computed in determinantal form, via the so-called Kasteleyn matrix K. This is a square matrix of size \(L^2\times L^2\) with rows/columns indexed by black/white vertices b/w of \({\mathbb {T}}_L\), as follows. If bw are not neighbors, then \(K(b,w)=0\). Otherwise, if (bw) is an edge of type r one sets \(K(b,w)=K_r\), cf. (2.28). We actually need four Kasteleyn matrices \(K_{{\varvec{\theta }}}\), \({\varvec{\theta }}=(\theta _1,\theta _2)\in \{0,1\}^2\), where the two indices label periodic/anti-periodic boundary conditions (depending on whether the index is 0/1) in the directions \({{e}}_i\). To obtain \(K_{{\varvec{\theta }}}\) from K, one multiplies by \((-1)^{\theta _1}\) (resp. by \((-1)^{\theta _2}\)) the matrix elements corresponding to edges (bw) where w has first coordinate equal L and b has first coordinate equal 1 (resp. w has second coordinate equal L and b has second coordinate equal 1). See Fig. 1. Of course, \(K_{00}=K\). We have then [32, 33] that

$$\begin{aligned} Z^0_L=\frac{1}{2}\sum _{{\varvec{\theta }}\in \{0,1\}^2}c_{{\varvec{\theta }}} \det (K_{{\varvec{\theta }}}) \end{aligned}$$
(3.1)

where \(c_{\varvec{\theta }}\in \{-1,+1\}\) and, moreover, three of the \(c_{\varvec{\theta }}\) have the same sign and the fourth one has the opposite sign. More precisely, for the square grid, with our choice of Kasteleyn matrix, one finds

(3.2)

(recall that we are assuming that L is even). The matrices \(K_{{\varvec{\theta }}}\) are diagonalized in the Fourier basis and

$$\begin{aligned} \det (K_{{\varvec{\theta }}})=\prod _{k\in {\mathcal {P}}({\varvec{\theta }})}\mu (k), \end{aligned}$$
(3.3)

where \(\mu (\cdot )\) is as in (2.5) and

$$\begin{aligned} {\mathcal {P}}({\varvec{\theta }})=\{k=(k_1,k_2),k_i=\frac{2\pi }{L} \left( n_i+\theta _i/2\right) ,n_i=0,\dots ,L-1\}. \end{aligned}$$
(3.4)

The matrices \(K_{{\varvec{\theta }}}\) are not necessarily invertible (e.g., if \(t_i\equiv 1\) then \(K_{00}\) is not because \(\mu (0)=0\)) and this question will play a role in Sect. 6. When the four matrices \(K_{{\varvec{\theta }}}\) are invertible, the correlation functions of the non-interacting measure can be written as

$$\begin{aligned}&{\mathbb {P}}_L(e_1,\dots ,e_k\in M)=\frac{1}{2Z^0_L}\nonumber \\&\quad \times \sum _{{\varvec{\theta }}\in \{0,1\}^2}c_{{\varvec{\theta }}}\det (K_{{\varvec{\theta }}}) \left[ \prod _{j=1}^k K_{{\varvec{\theta }}}(b_j,w_j)\right] \det \{K^{-1}_{{\varvec{\theta }}}(w_n,b_m)\}_{1\le n,m\le k} \end{aligned}$$
(3.5)

where the edge \(e_j\) has black/white vertex \(b_j/w_j\). The inverse of the matrix \(K_{{\varvec{\theta }}}\) can be computed explicitly as

$$\begin{aligned} K^{-1}_{{\varvec{\theta }}}(w_x,b_y)=\frac{1}{L^2}\sum _{k\in {\mathcal {P}}({\varvec{\theta }})}\frac{e^{-i k (x-y)}}{\mu (k)}=:g_{L}^{\varvec{\theta }}(x,y), \end{aligned}$$
(3.6)

where \(w_x\) (resp. \(b_y\)) is the white (resp. black) site with coordinate x (resp. y). Provided that

$$\begin{aligned} |k-p^\pm |\gg L^{-2}, \qquad \forall k\in {\mathcal {P}}({\varvec{\theta }}), \end{aligned}$$
(3.7)

it is easy to see that \(K^{-1}_{{\varvec{\theta }}}(w_x,b_y)=g(x,y)+o(1)\) as \(L\rightarrow \infty \), where

$$\begin{aligned} g(x,y):=\int _{[-\pi ,\pi ]^2}\frac{dk}{(2\pi )^2}\frac{e^{-i k(x-y)}}{\mu (k)}. \end{aligned}$$
(3.8)

Condition (3.7) can fail for some values of L and of \({\varvec{\theta }}\). For this reason, in Sect. 6 the values \(k^\pm _{\varvec{\theta }}\in {\mathcal {P}}({\varvec{\theta }})\) that are closest to the zeros of \(\mu \) will be treated separately, see in particular Sects. 6.1 and 6.5.

Due to the fact that \(\mu \) has two simple zeros, the matrix element g(xy) decays only as the inverse distance between \(w_x\) and \(b_y\). More precisely

$$\begin{aligned} g(x,y) = \frac{1}{2\pi }\sum _{\omega =\pm } \frac{e^{-ip^\omega (x-y)}}{\phi _\omega (x-y)}+r(x,y) \end{aligned}$$
(3.9)

where \(r(x,y)=O(1/|x-y|^2)\) and \(\phi _\omega \) was defined in (2.20).

3.2 Grassmann representation of the generating functions

We refer for instance to [21] for an introduction to Grassmann variables and Grassmann integration; here we just recall a few basic facts. It is well known that determinants can be represented as Gaussian Grassmann integrals. For our purposes, we associate a Grassmann variable \(\psi ^+_x\) (resp. \(\psi ^-_x\)) with the black (resp. white) site indexed x. We denote by \(\int D\psi f(\psi )\) the Grassmann integral of a function f and since the variables \(\psi ^\pm _x\) anti-commute among themselves and there is a finite number of them, we need to define the integral only for polynomials f. The Grassmann integration is a linear operation that is fully defined by the following conventions:

$$\begin{aligned} \int D\psi \, \prod _{x\in \Lambda }\psi ^-_x\psi ^+_x=1 , \end{aligned}$$
(3.10)

the sign of the integral changes whenever the positions of two variables are interchanged (in particular, the integral of a monomial where a variable appears twice is zero) and the integral is zero if any of the \(2|\Lambda |\) variables is missing. We also consider Grassmann intergrals of functions of the type \(f(\psi )=\exp (Q(\psi ))\), with Q a sum of monomials of even degree. By this, we simply mean that one replaces the exponential by its finite Taylor series containing only the terms where no Grassmann variable is repeated.

It is well known that the definition of Grassmann integration allows one to write the determinant of a matrix as the integral of the exponential of the associated Grassmann quadratic form (such integral will be called a “Gaussian Grassmann integral”, for the obvious formal analog with usual Gaussian integrals). In particular, recalling the definition of \(K_{{\varvec{\theta }}}\) given before (3.1),

$$\begin{aligned} \det (K_{{\varvec{\theta }}})=\int _{({\varvec{\theta }})} D\psi \, e^{S(\psi )}, \end{aligned}$$
(3.11)

where

$$\begin{aligned} S(\psi )=-\sum _{x,y\in \Lambda } K_{00}(b_x,w_y)\psi ^+_x\psi ^-_{y} \end{aligned}$$
(3.12)

and the index \(({\varvec{\theta }})\) below the integral means that one has to identify

$$\begin{aligned} \psi ^\pm _{(L+1,x_2)}:=(-1)^{\theta _1}\psi ^\pm _{(1,x_2)}, \quad \psi ^\pm _{(x_1,L+1)}:=(-1)^{\theta _2}\psi ^\pm _{(x_1,1)}. \end{aligned}$$
(3.13)

More compactly we write

$$\begin{aligned} S(\psi )=-\sum _e E_e \end{aligned}$$

where the sum runs over edges of \({\mathbb {T}}_L\) and, if e is an edge (bw),

$$\begin{aligned} E_e=K_{00}(b,w)\psi ^+_{x(b)}\psi ^-_{x(w)}. \end{aligned}$$
(3.14)

Our goal here is to express, via a Grassmann integral, the partition function of the interacting dimer model, and more generally the generating function \({\mathcal {W}}_\Lambda (A)\) defined by

$$\begin{aligned} e^{{\mathcal {W}}_\Lambda (A)}:=\sum _{M\in \Omega _{L}}p_{L,\lambda }(M)\prod _{e}e^{A_e {\mathbb {1}}_{e}{(M)}} \end{aligned}$$
(3.15)

where the product runs over the edges of \({\mathbb {T}}_L\) and \(A_e\in {\mathbb {R}}\). Note that \(e^{{\mathcal {W}}_L(0)}\) is the partition function and that any multi-dimer truncated correlation function of the type \({\mathbb {E}}_{L,\lambda } ( {\mathbb {1}}_{e_1};\dots ;{\mathbb {1}}_{e_k})\) can be obtained by differentiating \({\mathcal {W}}_\Lambda (A)\) with respect to \(A_{e_1},\dots ,A_{e_k}\) and setting \(A\equiv 0\).

Recall that the perturbed probability weight \(p_{L,\lambda }\) depends on the local ‘energy function’ f via (2.11)–(2.12). Without loss of generality, we can assume that (2.12) holds with

$$\begin{aligned} f(M)=\sum _{s=1}^n c_s{\mathbb {1}}_{P_s}(M) \end{aligned}$$
(3.16)

where \(c_s\) are real constants, n is an integer, \(P_s\) are finite collections of edges such that no space translation of \(P_s\) coincides with a \(P_{s'}, s\ne s'\) and \({\mathbb {1}}_{P_s}=\prod _{e\in P_s}{\mathbb {1}}_e\) is the indicator that all edges in \(P_s\) belong to M. Again without loss of generality we assume that each \(P_s\) contains at least 2 edges (if \(P_s\) consists in just one edge, its effect is just to modify the weights \({\underline{t}}\)). Under these assumptions, the following representation holds.

Proposition 1

Let \(\lambda \) be small enough. Then, one has

$$\begin{aligned} e^{{\mathcal {W}}_L(A)}= \frac{1}{2}\sum _{{\varvec{\theta }}\in \{0,1\}^2}c_{{\varvec{\theta }}} \int _{({\varvec{\theta }})} D\psi \,e^{S(\psi )+V(\psi ,A)} \end{aligned}$$
(3.17)

where

$$\begin{aligned} V(\psi ,A)=-\sum _e (e^{A_e}-1)E_e+\sum _{\gamma \subset \Lambda }c(\gamma )\prod _{b\in \gamma }E_b e^{A_b}. \end{aligned}$$
(3.18)

The first sum runs over all edges of \({\mathbb {T}}_L\) and \(E_e\) is as in (3.14). In the second sum, \(\gamma \) are finite subsets of disjoint edges of \({\mathbb {T}}_L\) such that \(|\gamma |\ge 2\), and \(c(\gamma )\) is a real constant satisfying translation invariance (\(c(\gamma )=c(\tau _x\gamma )\)) and the bound

$$\begin{aligned} |c(\gamma )|\le (a |\lambda |)^{\max \{1,b \delta (\gamma )\}}, \end{aligned}$$
(3.19)

for some constants \(a,b>0\), independent of L, and \(\delta (\gamma )\) the tree distance of \(\gamma \), that is, the length of the shortest tree graph on \(\Lambda \) containing \(\gamma \) (the precise definition of \(c(\gamma )\) is given below).

Remark 6

Both \(S(\psi )\) and \(V(\psi ,A)\) are invariant under the following symmetry transformation of the Grassmann fields:

$$\begin{aligned} \psi ^\pm _x\rightarrow (-1)^x\psi ^\pm _x,\qquad c\rightarrow c^*, \end{aligned}$$
(3.20)

where \(c\rightarrow c^*\) indicates that all the constants appearing in \(S(\psi )\) and \(V(\psi ,A)\) are mapped to their complex conjugates. Also, we used the notation \((-1)^x:=(-1)^{x_1+x_2}\). It is straightforward to check that, under this transformation, \(E_e\rightarrow E_e\), for all the edges e, which clearly shows that the considered transformation is in fact a symmetry of the Grassmann action. This symmetry will play a role in Sect. 6, in reducing the number of independent running coupling constants arising in the multiscale computation of the Grassmann generating function.

Proof of Proposition 1

The proposition has been proven in [25] in the case of constant weights \(t_i\equiv 1\) and plaquette interaction as in (2.14); the extension to the present situation is rather straightforward, so we will be concise.

Let

$$\begin{aligned} S=\{\tau _x P_s,s=1,\dots ,n, x\in \Lambda \} \end{aligned}$$

and remark that by assumption all elements of S are distinct and contain at least two edges. If \(B\in S\) is a space translation of \(P_s\), set

$$\begin{aligned} u(B)=\exp (\lambda c_s)-1. \end{aligned}$$
(3.21)

We start by writing

$$\begin{aligned} e^{{\mathcal {W}}_L(A)}= & {} \sum _{M\in \Omega _L} w^{(A)}(M)\prod _{x\in \Lambda }\prod _{s=1}^n(1+(e^{\lambda c_s}-1) {\mathbb {1}}_{\tau _x P_s}(M))\nonumber \\= & {} Z_L^{0,(A)} \sum _{\sigma \subset S}{\mathbb {E}}_{L}^{(A)}\left( \prod _{B\in \sigma }u(B){\mathbb {1}}_B(M)\right) \end{aligned}$$
(3.22)

with

$$\begin{aligned} w^{(A)}(M)=t_1^{N_1(M)}t_2^{N_2(M)}t_3^{N_3(M)}e^{\sum _{b\in M}A_b}, \qquad Z_L^{0,(A)}=\sum _{M\in \Omega _L} w^{(A)}(M) \end{aligned}$$

and \({\mathbb {P}}_L^{(A)}\) the probability measure with density \(w^{(A)}(M)/Z_L^{0,(A)}\). By manipulating the sum in the r.h.s. of (3.22), one can rewrite it as

$$\begin{aligned} \sum _{n\ge 0}\sum ^*_{\gamma _1,\dots ,\gamma _n} Z_L^{0,(A)}{\mathbb {E}}_{L}^{(A)}\left( \prod _{i=1}^n {{\tilde{c}}}(\gamma _i){\mathbb {1}}_{\gamma _i}(M)\right) \end{aligned}$$
(3.23)

where the term \(n=0\) has to be interpreted as equal to 1 and the sum \(\sum ^*\) is over non-empty, mutually disjoint subsets \(\gamma _i\) of edges of \({\mathbb {T}}_L\). The constant \({{\tilde{c}}}(\gamma )\) is given as follows. Let \(\Sigma _\gamma \) be the set of all collections of the type \(Y=\{B_1,\dots ,B_{|Y|}\}\) where: \(B_i\in S\), \(B_i\ne B_j\) for \(i\ne j\), \(\cup _i B_i=\gamma \) and such that Y cannot be divided into two non-empty sub-collections \(\{B_{i_1},\dots ,B_{i_k}\}\) and \(\{B_{i_{k+1}},\dots ,B_{i_{|Y|}}\}\) with \((\cup _{j\le k}B_{i_j})\cap (\cup _{j>k} B_{i_j})=\emptyset \). Then

$$\begin{aligned} {{\tilde{c}}}(\gamma )=\sum _{Y\in \Sigma _\gamma }\prod _{B\in Y} u(B). \end{aligned}$$
(3.24)

Now we rewrite (3.23) as

$$\begin{aligned} \sum _{n\ge 0}\sum ^*_{\gamma _1,\dots ,\gamma _n} \prod _{j=1}^n {{\tilde{c}}}(\gamma _j) \left[ \prod _{b\in \gamma _j}\partial _{A_b}\right] Z_L^{0,(A)}. \end{aligned}$$
(3.25)

The partition function \( Z_L^{(A)}\) corresponds to a non-interacting dimer model, with edge-dependent weights \(t_e e^{A_e}\). Then, as in (3.1) and (3.11) we have

$$\begin{aligned} Z_L^{0,(A)}=\frac{1}{2}\sum _{{\varvec{\theta }}\in \{0,1\}^2}c_{{\varvec{\theta }}} \int _{({\varvec{\theta }})}D\psi \, e^{S(\psi )-\sum _e (e^{A_e}-1)E_e}. \end{aligned}$$
(3.26)

Using expression (3.26) in (3.25) one readily concludes, as in [25], that (3.18) holds with

$$\begin{aligned} c(\gamma )=(-1)^{|\gamma |}{{\tilde{c}}}(\gamma ). \end{aligned}$$
(3.27)

If \(\lambda \) is small enough, it is easy to see that the bound (3.19) holds. \(\square \)

For the 6-vertex model with interaction (2.15), the potential V is exactly quartic in the fields \(\psi \): indeed, \(c(\gamma )\ne 0\) only if \(\gamma \) is the pair of edges \(\gamma =\{e_1,e_{2}\}\) or \(\gamma =\{e_{3},e_{4}\}\) as in Fig. 2 or a translation thereof. For the plaquette model with interaction (2.14), instead, \(c(\gamma )\) is non-zero only if \(\gamma \) is a collection of \(|\gamma |\ge 2\) adjacent parallel edges, in which case \(c(\gamma )=(-1)^{|\gamma |}(e^\lambda -1)^{|\gamma |-1}\).

In the following (in the comparison between the discrete lattice model and the continuum reference model) we will also need the generating function for mixed dimer and fermionic correlations. Namely, let \(\{\phi ^+_x,\phi ^-_x\}_{x\in \Lambda }\) be Grassmann variables that anti-commute among themselves and with the \(\psi ^\pm \) variables. Then, we let

$$\begin{aligned} e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,\phi )}:=\int _{({\varvec{\theta }})} D\psi \, e^{S(\psi )+V(\psi ,A)+(\psi ,\phi )} \end{aligned}$$
(3.28)

and

$$\begin{aligned} e^{{\mathcal {W}}_L(A,\phi )}:= \frac{1}{2}\sum _{{\varvec{\theta }}\in \{0,1\}^2}c_{{\varvec{\theta }}} e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,\phi )}. \end{aligned}$$
(3.29)

Here, \(V(\psi ,A)\) is as in Proposition 1, while

$$\begin{aligned} (\psi ,\phi ):=\sum _{x\in \Lambda } (\psi ^+_x\phi ^-_x+\phi ^+_x\psi ^-_x). \end{aligned}$$

We define \(g_{L}(e_1,\dots ,e_k;x_1,\dots ,x_n;y_1,\dots ,y_n)\) as the truncated correlations associated with the generating functionFootnote 5\({\mathcal {W}}_L(A,\phi )\):

$$\begin{aligned}&g_{L}(e_1,\dots ,e_k;x_1,\dots ,x_n;y_1,\dots ,y_n)\nonumber \\&\quad :=\left. \partial _{A_{e_1}}\dots \partial _{A_{e_k}}\partial _{\phi ^-_{y_1}}\dots \partial _{\phi ^-_{y_n}} \partial _{\phi ^+_{x_1}}\dots \partial _{\phi ^+_{x_n}} {\mathcal {W}}_{L}(A,\phi )\right| _{A\equiv 0,\phi \equiv 0}. \end{aligned}$$
(3.30)

Two cases that will play a central role in the following are \(k=0\), \(n=1\) (the interacting propagator), and \(k=n=1\) (the interacting vertex function), which deserve a distinguished notation.

Interacting propagator:

$$\begin{aligned} g_{L}(\emptyset ;x;y)=\frac{1}{2Z_{L}}\sum _{{\varvec{\theta }}}c_{{\varvec{\theta }}}\int _{({\varvec{\theta }})}D \psi \, e^{S(\psi )+V(\psi ,0)}\psi ^-_{x}\psi ^+_{y}=:G^{(2)}_L(x,y); \end{aligned}$$
(3.31)

that is, \(G^{(2)}_L(x,y)={\langle \psi ^-_{x}\psi ^+_{y}\rangle }_L\), where \({\langle f\rangle }_L\) indicate the Grassmann “average” \(\frac{1}{2Z_{L}}\sum _{{\varvec{\theta }}}c_{{\varvec{\theta }}}\int _{({\varvec{\theta }})}D\psi \, e^{S(\psi )+V(\psi ,0)}f(\psi )\).

Interacting vertex function:

if \({\mathcal {I}}_{e}=\partial _{A_e}V(\psi ,A)\big |_{A=0}\) is the Grassmann counterpart of the dimer observable at e, and e is an edge of type r with black site labelled z, then

$$\begin{aligned} g_{L}(e;x;y)={\langle {{\mathcal {I}}}_{e}\psi ^-_{x}\psi ^ +_{y}\rangle }_L-{\langle {{\mathcal {I}}}_{e}\rangle }_L{\langle \psi ^-_{x}\psi ^+_{y}\rangle }_L =:G^{(2,1)}_{r,L}(z,x,y); \end{aligned}$$
(3.32)

that is, \(G^{(2,1)}_{r,L}(z,x,y)={\langle {\mathcal {I}}_e;\psi ^-_x\psi ^+_y\rangle }_L\), where the semicolon indicates truncated expectation.

In the following we will also need a distinguished notation for the two-point dimer–dimer correlation: if \(e_1,e_2\) are two edges of type \(r,r'\), and black sites labelled xy, respectively, we let

$$\begin{aligned} g_{L}(e_1,e_2;\emptyset ;\emptyset )=:G^{(0,2)}_{r,r',L}(x,y). \end{aligned}$$
(3.33)

Note that all the correlations \(g_L(e_1,\dots ,e_k;x_1,\dots ,x_n;y_1,\dots ,y_n)\) are well defined for any finite L, despite the fact that the Kasteleyn matrix \(K_{{\varvec{\theta }}}\) may not be invertible for some choices of \({\varvec{\theta }},L\). The multipoint correlations,

$$\begin{aligned} g_{L}(e_1,\dots ,e_k;x_1,\dots ,x_n;y_1,\dots ,y_n), \end{aligned}$$

admit a thermodynamic limit as \(L\rightarrow \infty \), as shown in Sect. 6; the limit can be expressed as a convergent multiscale fermionic expansion and will be denoted

$$\begin{aligned} g(e_1,\ldots ,e_k;x_1,\dots ,x_n;y_1,\dots ,y_n). \end{aligned}$$

In particular, the thermodynamic limit of the two-point dimer–dimer correlation will be denoted by \(G^{(0,2)}_{r,r'}(x,y)\), while the \(L\rightarrow \infty \) limit of the interacting propagator and vertex function will be denoted \(G^{(2)}(x,y)\) and \(G^{(2,1)}_r(z,x,y)\).

3.3 Lattice Ward Identity

The generating function \({\mathcal {W}}_{L}(A,\phi )\) has a gauge symmetry property that implies certain identities (lattice Ward identities) involving its derivatives. These identities were derived in [26] for the model with \(t_i\equiv 1\) and they hold (with the same proof) also for the general model studied here. We recall here, without giving the proof, the Ward Identity for the ‘vertex function’, but similar relations can be easily derived for higher point correlations: for any finite L,

$$\begin{aligned}&\sum _{r=1}^4 G^{(2,1)}_{r,L}(x,y,z)=-\delta _{x,z}G^{(2)}_L(y,x), \end{aligned}$$
(3.34)
$$\begin{aligned}&\sum _{r=1}^4 G^{(2,1)}_{r,L}(x-v_r,y,z)=-\delta _{x,y}G^{(2)}_L(x,z), \end{aligned}$$
(3.35)

with \(\delta _{x,y}\) the Krokecker delta, see [26, Eq.(4.9)–(4.10)]. By taking the difference between these two equations, we get (see [26, Eq.(4.17)])

$$\begin{aligned} \delta _{x,y}G^{(2)}_L(x,z) -\delta _{x,z}G^{(2)}_L(y,x)=-\sum _{r=2}^4\nabla _{-v_r}G^{(2,1)}_{r,L}(x,y,z), \end{aligned}$$
(3.36)

where \((\nabla _n f)(x,y,z):=f(x+n,y,z)-f(x,y,z)\) is the (un-normalized) discrete derivative acting on the x variable. By taking the limit \(L\rightarrow \infty \), we see that (3.34)–(3.36) also hold for the infinite volume correlation functions \(G^{(2)}(x,y), G^{(2,1)}_{r}(x,y,z)\).

In Fourier space, we define

$$\begin{aligned} {{\hat{G}}}^{(2)}(p)= & {} \sum _{x} G^{(2)}(x,0)e^{i p x} \end{aligned}$$
(3.37)
$$\begin{aligned} {{\hat{G}}}^{(2,1)}_r(k,p)= & {} \sum _{x,z}e^{-i px-ikz} G^{(2,1)}_r(x,0,z) \end{aligned}$$
(3.38)
$$\begin{aligned} {{\hat{G}}}^{(0,2)}_{r,r'}(p)= & {} \sum _x e^{-ipx} G^{(0,2)}_{r,r'}(x,0). \end{aligned}$$
(3.39)

Then, the infinite-volume limit of (3.34)–(3.36) can be rewritten as

$$\begin{aligned}&\sum _{r=1}^4 {{\hat{G}}}^{(2,1)}_{r}(k,p)=-{{\hat{G}}}^{(2)}(k+p), \end{aligned}$$
(3.40)
$$\begin{aligned}&{{\hat{G}}}^{(2)}(k+p)-{{\hat{G}}}^{(2)}(k)=\sum _{r=2}^4(e^{-ipv_r}-1) {{\hat{G}}}^{(2,1)}_r(k,p). \end{aligned}$$
(3.41)

In the following the asymptotic behavior at large distances of the interacting propagator and vertex function will be computed in terms of a reference continuum model, see next section, which plays the role of the ‘infrared fixed point’ of our lattice dimer model in its Grassmann formulation.

4 The Infrared Fixed Point Theory

In order to introduce the “infra-red fixed point” of our theory (referred to in the following as “the continuum model” or “the reference model”), we need a couple of preliminary definitions. First, we let \({\mathcal {M}}\) be the \(2\times 2\) matrix with unit determinant

$$\begin{aligned} {\mathcal {M}}=\frac{1}{\sqrt{\Delta }} \begin{pmatrix} {\bar{\beta }}^1&{} {\bar{\beta }}^2\\ -{\bar{\alpha }}^1&{}-{\bar{\alpha }}^2 \end{pmatrix} \end{aligned}$$
(4.1)

where \({\bar{\alpha }}^j,{\bar{\beta }}^j\in {\mathbb {R}}\), \(j=1,2\) and \(\Delta ={\bar{\alpha }}^1{\bar{\beta }}^2-{\bar{\alpha }}^2{\bar{\beta }}^1>0\) (for the moment, these are free parameters; eventually, they will be the real and imaginary parts of the functions \({\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \) that appear in Theorem 1). Also, given \(L>0\) (the system size), an integer N (ultra-violet cut-off) and \(Z>0\), we introduce a Grassmann Gaussian integrationFootnote 6\(P_Z^{[\le N]}(d\psi )\) on the family of Grassmann variables

$$\begin{aligned} \{ {\hat{\psi }}^\pm _{k,\omega },\; \omega =\pm 1, \; k\in {\mathcal {K}}\},\qquad {\mathcal {K}}=\left\{ {\mathcal {M}}\cdot p\; \left| \;p\in \left( \frac{2\pi }{L}\right) ({\mathbb {Z}}+1/2)^2\right. \right\} , \end{aligned}$$

defined by the propagator

$$\begin{aligned} \int P_Z^{[\le N]}(d\psi ){\hat{\psi }}^-_{k,\omega }\psi ^+_{k',\omega '}= \delta _{\omega ,\omega '}\delta _{k,k'}\frac{L^2 }{Z} \frac{\chi _{N}(k)}{ {\bar{D}}_\omega (k)} \end{aligned}$$
(4.4)

where:

  • \(\chi _{N}(k)=\chi (2^{-N}|{\mathcal {M}}^{-1} k|)\), with \(\chi :{\mathbb {R}}^+\rightarrow [0,1]\) a \(C^\infty \) cut-off function that is equal to 1 if its argument is smaller than 1 and equal to 0 if its argument is larger than 2;

  • \({\bar{D}}_\omega (k)={\bar{\alpha }}_\omega k_1+{\bar{\beta }}_\omega k_2\), with

    $$\begin{aligned} {\bar{\alpha }}_\omega =\omega {\bar{\alpha }}^1+i {\bar{\alpha }}^2,\quad {\bar{\beta }}_\omega =\omega {\bar{\beta }}^1+i {\bar{\beta }}^2. \end{aligned}$$
    (4.5)

    Observe that, since we are assuming \(\Delta >0\), we have that

    $$\begin{aligned} \frac{{\bar{\alpha }}_\omega }{{\bar{\beta }}_\omega }\not \in {\mathbb {R}}. \end{aligned}$$
    (4.6)

While \({\mathcal {K}}\) is an infinite set, we effectively have only a finite number of non-zero Grassmann variables \({\hat{\psi }}^\pm _{k,\omega }\), because \(\chi _N(k)\) is non-zero only for a finite number of values of k in \({\mathcal {K}}\).

Note that, setting \(q={\mathcal {M}}^{-1}k\), the r.h.s. of (4.3) equals

$$\begin{aligned} \delta _{\omega ,\omega '}\delta _{q,q'}\frac{L^2}{Z \sqrt{\Delta }}\frac{\chi (2^{-N}|q|)}{-i q_1+\omega q_2}. \end{aligned}$$
(4.7)

In the language of Quantum Field Theory, in the limit \(\lim _{L\rightarrow \infty ,N\rightarrow \infty }\), (4.6) is just the propagator of chiral massless relativistic fermions.

It is convenient to define, for \(x\in {\mathbb {R}}^2\), the Grassmann variables

$$\begin{aligned} \psi ^\pm _{x,\omega }:=\frac{1}{L^2}\sum _{k\in {\mathcal {K}}}e^{\pm i k x}{\hat{\psi }}^\pm _{k,\omega }. \end{aligned}$$
(4.2)

Note that \(\psi ^\pm _{x,\omega }\) has anti-periodic boundary conditions on

$$\begin{aligned} \Lambda :=({\mathcal {M}}^T)^{-1}{\mathcal {T}}_L, \quad {\mathcal {T}}_L={\mathbb {R}}^2/(L{\mathbb {Z}}^2) \end{aligned}$$

and that

$$\begin{aligned} \frac{g^{[\le N]}_{R,\omega }(x-y)}{Z}:=\int P_Z^{[\le N]}(d\psi )\psi ^-_{x,\omega }\psi ^+_{y,\omega }=\frac{1}{ ZL^2}\sum _{k\in {\mathcal {K}}}e^{-i k(x-y)}\frac{\chi _{N}(k)}{{\bar{D}}_\omega (k)}. \end{aligned}$$
(4.8)

The generating functional \({\mathcal {W}}_{L,N}(J,\phi )\) of the continuum model is

$$\begin{aligned} e^{{{\mathcal {W}}}_{L,N}(J,\phi )} = {\int \! P_Z^{[\le N]}(d\psi ) e^{{{\mathcal {V}}}(\sqrt{Z}\psi ) + \sum _{j=1}^2 (J^{(j)},\,\rho ^{(j)})+ Z\,(\psi ,\phi )}} \;, \end{aligned}$$
(4.9)

where \(J=\{J^{(j)}_{x,\omega }\}^{j=1,2}_{\omega =\pm ,\, x\in \Lambda }\) are external “sources” (real-valued test functions) and \(\phi =\{\phi ^\sigma _{x,\omega }\}^{\sigma ,\omega =\pm }_{x\in \Lambda }\) are “external Grassmann sources”, i.e. \(\phi ^\sigma _{x,\omega }\) is a Grassmann variable. Also, we used the notation

$$\begin{aligned} (J^{(j)},\rho ^{(j)}):=\sum _{\omega =\pm }\int _\Lambda dx\; J^{(j)}_{x,\omega }\rho ^{(j)}_{x,\omega }, \end{aligned}$$

with

$$\begin{aligned} \rho ^{(1)}_{x,\omega } = \psi ^+_{x,\omega } \psi ^{-}_{x,\omega }\;, \qquad \rho ^{(2)}_{x,\omega } = \psi ^+_{x,\omega } \psi ^{-}_{x,-\omega }\; \end{aligned}$$
(4.10)

and

$$\begin{aligned} (\psi ,\phi ):=\sum _{\omega =\pm }\int _\Lambda dx\,(\psi ^+_{x,\omega }\phi ^-_{x,\omega }+\phi ^+_{x,\omega }\psi ^-_{x,\omega })\;. \end{aligned}$$

Finally, the interaction \({\mathcal {V}}\) in (4.9) is

$$\begin{aligned} {{\mathcal {V}}}(\psi )=\frac{\lambda _\infty }{2} \sum _{\omega =\pm }\int _\Lambda dx\int _\Lambda dy\; v(x-y) \psi ^+_{x,\omega } \psi ^-_{x,\omega }\psi ^+_{y,-\omega }\psi ^-_{y,-\omega }\;, \end{aligned}$$
(4.11)

where \(\lambda _\infty \in {\mathbb {R}}\), \(v(x)=v_0({\mathcal {M}}^T x)\) and \(v_0(\cdot )\) is a smooth rotationally invariant potential, exponentially decaying to zero at large distances, normalized as

$$\begin{aligned} \int _{{\mathbb {R}}^2}dx\, v_0(x)= \int _{{\mathbb {R}}^2}dx\, v(x)=1. \end{aligned}$$
(4.12)

We emphasize that, while (4.11) seems to depend on an uncountable set of Grassmann variables \(\{\psi ^\pm _{x,\omega },\phi ^\pm _{x,\omega }\}_{x\in \Lambda }\), writing everything in Fourier space there is only a finite number of non-zero Grassmann variables.

In the special case \({\bar{\alpha }}_\omega =(-i-\omega ),{\bar{\beta }}_\omega =(-i+\omega )\), which is relevant for the interacting dimer model with \({\underline{t}}\equiv 1\), the continuum model reduces to that studied in [26, Sec. 5], if the constants \(Z^{(1)}\) and \(Z^{(2)}\) that appear there are fixed to 1. Setting instead \({\bar{\alpha }}_\omega =-i, {\bar{\beta }}_\omega =\omega \) in (4.9) (so that \(\Delta =1\)) one obtains, apart from minor differences, the model studied in [6, Sec. 3] and [12, Sec. 3].

Remark 7

In order to recognize the equivalence of the model (4.9) with \({\bar{\alpha }}_\omega =-i, {\bar{\beta }}_\omega =\omega \) and the one in, e.g., [12, Section 3] (or, analogously, the one in [6, Section 3]), one needs to set to zero some of the external fields, rotate the coordinate system and rescale some constants. More precisely, if \({\mathbb {W}}_{L,N}(J^{(1)},\phi )\) denotes the generating functional used in [12] with \(J^{(1)}_{x,\omega }=Z^{(3)}J_x+\omega {{\tilde{Z}}}^{(3)}{{\tilde{J}}}_x\), see [12, Eq. (28)], then, setting \(J^{(2)}_x\equiv 0\) in (4.9),

$$\begin{aligned} {\mathcal {W}}_{L,N}((J^{(1)},0),\phi ;\lambda _\infty ) =const.+{\mathbb {W}}_{L,N}({\mathcal {J}}^{(1)},\varphi ; -\Delta ^{-1}\lambda _\infty ) \end{aligned}$$
(4.13)

where the constant is independent of \(J^{(1)},\phi \) (so that it does not influence the correlation functions; it depends upon \(\Delta \) and is due to the rescaling of the Grassmann fields), while

$$\begin{aligned} {\mathcal {J}}^{(1)}(x):=\Delta ^{1/2}J^{(j)}(({\mathcal {M}}^T)^{-1}x), \quad \varphi ^\pm (x):=\Delta ^{1/4}\phi ^\pm (({\mathcal {M}}^T)^{-1}x), \end{aligned}$$
(4.14)

and we denoted explicitly the dependence of the generating function on \(\lambda _\infty \). This immediately implies obvious relations between the correlation functions \(G^{(2,1)}_{R,\omega ',\omega }(x,y,z)\), \(G^{(2)}_{R,\omega }(x,y)\) and \(S^{(j,j)}_{R,\omega ,\omega '}(x,y)\), defined below, and the analogous ones of [12].

The peculiarity of the continuum model is that its correlations can be computed exactly, in the removed cutoff limit where \(N\rightarrow \infty \) first and then \(L\rightarrow \infty \). This is because, as compared to its lattice counterpart, the continuum model is “chiral gauge invariant”, which means that the correlation functions satisfy two hierarchies of Ward Identities, distinguished by the choice of the ‘chirality index’ \(\omega \), see (4.24) below. These additional symmetries, together with other identities among correlation functions (the so-called Schwinger-Dyson equations), allow one to get closed equations for correlation functions, in the removed cutoff limit. In this sense, the infrared fixed point theory can be regarded as “integrable”.

We define the following correlation functions of the reference model: if xyz are distinct points of \(\Lambda \),

$$\begin{aligned}&G^{(2,1;L,N)}_{R,\omega ',\omega }(x,y,z) = \frac{\partial ^3}{\partial J_{x,\omega '}^{(1)}\partial \phi ^-_{z,\omega }\partial \phi ^+_{y,\omega } } {{\mathcal {W}}}_{L,N}(J,\phi )|_{J=\phi =0} \nonumber \\&G^{(2;L,N)}_{R,\omega }(x,y):=\frac{\partial ^2}{ \partial \phi ^-_{y,\omega } \partial \phi ^+_{x,\omega }} {{\mathcal {W}}}_{L,N}(J,\phi )|_{J=\phi =0} \nonumber \\&S^{(j,j;L,N)}_{R,\omega ,\omega '}(x,y):=\frac{\partial ^2}{\partial J_{x,\omega }^{(j)} \partial J_{y,\omega '}^{(j)}}{{\mathcal {W}}}_{L,N}(J,\phi )|_{J=\phi =0}. \end{aligned}$$
(4.15)

From the construction of the correlation functions of the model, see e.g. [6, Section 3 and 4], one obtains in particular the existence of the following limits where cut-offs are removed:

$$\begin{aligned}&G^{(2,1)}_{R,\omega ',\omega }(x,y,z) =\lim _{L\rightarrow \infty }\lim _{N\rightarrow \infty } G^{(2,1;L,N)}_{R,\omega ',\omega }(x,y,z)\;,\nonumber \\&G^{(2)}_{R,\omega }(x,y) =\lim _{L\rightarrow \infty }\lim _{N\rightarrow \infty }G^{(2;L,N)}_{R,\omega }(x,y) \;,\nonumber \\&S^{(j,j)}_{R,\omega ,\omega '}(x,y)=\lim _{L\rightarrow \infty } \lim _{N\rightarrow \infty }S^{(j,j;L,N)}_{R,\omega ,\omega '}(x,y)\;. \end{aligned}$$
(4.16)

Away from \(x=0\), the so-called “density-density” correlation \(S^{(1,1)}\) is given by [26, Eq. (5.12)]

$$\begin{aligned} S^{(1,1)}_{R,\omega ,\omega }(x,0) = \frac{1}{4\pi ^2 Z^2(1-\tau ^2)}\frac{1}{({\bar{\phi }}_\omega (x))^2}+R_1(x), \end{aligned}$$
(4.17)

where \( {\bar{\phi }}_\omega (x):=\omega ({\bar{\beta }}_\omega x_1-{\bar{\alpha }}_\omega x_2)\),

$$\begin{aligned} \tau =-\frac{\lambda _\infty }{4\Delta \pi } \end{aligned}$$
(4.18)

and \(|R_1(x)|\le C|x|^{-3}\). On the other hand, the “mass-mass correlation” \(S^{(2,2)}\) satisfies (see [26, Eq.(6.14)])

$$\begin{aligned} S^{(2,2)}_{R,\omega ,-\omega }(x,0) = \frac{{\bar{B}}}{4\pi ^2Z^2}\frac{1}{|{\bar{\phi }}_\omega (x)|^{2\nu }}+R_2(x), \end{aligned}$$
(4.19)

where \({\bar{B}}\) is an analytic function of \(\lambda _\infty , Z, {\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \), which is equal to 1 at \(\lambda _\infty =0\),

$$\begin{aligned} \nu =\frac{1-\tau }{1+\tau }, \end{aligned}$$
(4.20)

see [26, Eq.(6.15)] and [26, Appendix C], and \(R_2\) is a correction term such that \(|R_2(x)|\le C|x|^{-2-\theta }\), for some \(\theta >0\) that, e.g., can be chosen \(\theta =1/2\).

We will not need the explicit form of \(G^{(2)}_{R,\omega }(x,0)\) and \( G^{(2,1)}_{R,\omega ',\omega }(x,0,z)\); let us just mention that they diverge as xz tend to zero but they are locally integrable functions (see, e.g., the expression of the interacting propagator in [6, eq.(4.18)]) and therefore admit Fourier transforms in the sense of distributionsFootnote 7,

$$\begin{aligned}&{{\hat{G}}}^{(2,1)}_{R,\omega ',\omega }(k,p)= \int dx \int dy\, e^{-ipx +i(k+p)y}\, G^{(2)}_{R,\omega ',\omega }(x,y,0) \;,\nonumber \\&{{\hat{G}}}^{(2)}_{R,\omega }(k) = \int dx\, e^{ikx}\,G^{(2)}_{R,\omega }(x,0). \end{aligned}$$
(4.21)

For later use, let us mention that the small-momenta behavior of \({{\hat{G}}}^{(2,1)}_{R,\omega ,\omega }\) and \({{\hat{G}}}^{(2)}_{R,\omega }\) are (see, e.g., [10, Theorem 2])

$$\begin{aligned} |{{\hat{G}}}^{(2)}_{R,\omega }(p)|\sim const\times {|p|}^{-1+O(\lambda _\infty ^2)}, \quad {\text {as}\quad p\rightarrow 0}, \end{aligned}$$
(4.22)

and, if \(0<{\mathfrak {c}}\le |p|,|k|,|k+p|\le 2{\mathfrak {c}}\),

$$\begin{aligned} |{{\hat{G}}}^{(2,1)}_{R,\omega ,\omega }(k,p)|\sim const\times {\mathfrak {c}}^{-2+O(\lambda _\infty ^2)}, \quad {\text {as}\quad {\mathfrak {c}}\rightarrow 0}. \end{aligned}$$
(4.23)

A very useful consequence of the exact solution of the continuum model is that the “propagator” and the “vertex function” satisfy the following Ward Identity (see [26, Eq.(5.9)]):

$$\begin{aligned} {Z}\sum _{\omega '=\pm }{\bar{D}}_{\omega '}(p){{\hat{G}}}^{(2,1)}_{R,\omega ',\omega }(k,p)= \frac{1}{1-\tau {{\hat{v}}}(p)} [{{\hat{G}}}^{(2)}_{R,\omega }(k) - {{\hat{G}}}^{(2)}_{R,\omega }(k+p)]\;. \end{aligned}$$
(4.24)

Note that this identity resembles formally the lattice Ward identity (3.41) of the dimer model, with the crucial difference that (4.24) are actually two identities (one for each choice of \(\omega \)).

5 Comparison Between Lattice and Continuum Model, and Proof of Theorems 12

The reason why the continuum model plays the role of the “infrared fixed point theory” for our interacting dimer model is that the large distance behaviour of the dimer correlation functions can be expressed in terms of linear combinations of the correlations of the continuum model, for a suitable choice of the parameters \(Z,\lambda _\infty ,{\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \). Let us spell out the explicit relation between correlation functions of the two models, in the special cases of the dimer interacting propagator, the vertex function and the dimer–dimer correlation. The result is a consequence of the multi-scale analysis described in Sect. 6 where, in particular, we prove the following:

Proposition 2

In a small neighborhood of \(\lambda =0\) in the complex plane, there exist:

  1. (1)

    two real-valued analytic functionsFootnote 8\(\lambda \mapsto {\bar{p}}^\omega \), with \(\omega =\pm \), called the interacting Fermi points, satisfying (2.42) and (2.45),

  2. (2)

    four complex-valued analytic functions \(\lambda \mapsto {\bar{\alpha }}_\omega ,\lambda \mapsto {\bar{\beta }}_\omega \) satisfying (2.41), (2.43),

  3. (3)

    two real-valued analytic functions \(\lambda \mapsto Z,\lambda \mapsto \lambda _\infty \) satisfying \(Z= 1+O(\lambda )\) and \(\lambda _\infty =O(\lambda )\),

such that the dimer–dimer correlation can be represented in the following form:

$$\begin{aligned}&G^{(0,2)}_{r,r'}(x,y) = \sum _{\omega =\pm } {{\hat{K}}}_{\omega ,r} {{\hat{K}}}_{\omega ,r'} S^{(1,1)}_{R,\omega ,\omega }(x,y)\nonumber \\&\qquad + \sum _{\omega =\pm } e^{i({\bar{p}}^\omega -{\bar{p}}^{-\omega })(x-y)} {\hat{H}}_{-\omega ,r} {{\hat{H}}}_{\omega ,r'} S^{(2,2)}_{R,\omega ,-\omega }(x,y)+R_{r,r'}(x,y)\;, \end{aligned}$$
(5.1)

where: \(\lambda \mapsto {{\hat{K}}}_{\omega ,r}\) and \(\lambda \mapsto {{\hat{H}}}_{\omega ,r}\) are complex-valued analytic function of \(\lambda \) satisfying \({{\hat{K}}}_{+,r}={{\hat{K}}}^*_{-,r}\), \({{\hat{H}}}_{+,r}={{\hat{H}}}^*_{-,r}\), \({{\hat{K}}}_{\omega ,r}=K_{\omega ,r}+O(\lambda )\), and \({{\hat{H}}}_{\omega ,r}=K_{\omega ,r}+O(\lambda )\); the correction term \(R_{r,r'}(x,y)\) is translational invariant and satisfies \(|R_{r,r'}(x,0)|\le C |x|^{-5/2}\).

Moreover, the interacting dimer propagator satisfies

$$\begin{aligned} {{\hat{G}}}^{(2)}(k+{\bar{p}}^\omega ) {\mathop {=}\limits ^{{\mathfrak {c}}\rightarrow 0}} {{\hat{G}}}^{(2)}_{R,\omega }(k)[1+O(|k|^\theta )], \end{aligned}$$
(5.2)

for some \(\theta >0\)Footnote 9. Finally, if \(0<{\mathfrak {c}}\le |p|,|k|,|k+p|\le 2{\mathfrak {c}}\), then the interacting vertex function of the dimer model satisfies

$$\begin{aligned} {{\hat{G}}}^{(2,1)}_r(k+ {\bar{p}}^\omega , p) {\mathop {=}\limits ^{{\mathfrak {c}}\rightarrow 0}} -\sum _{\omega '=\pm }{{\hat{K}}}_{\omega ',r}{{\hat{G}}}^{(2,1)}_{R,\omega ',\omega }(k,p) [1+O({\mathfrak {c}}^\theta )]\;. \end{aligned}$$
(5.3)

The statements contained in this proposition are proved in different subsections of Sect. 6: the construction of the functions \({\bar{p}}^\omega , {\bar{\alpha }}_\omega , {\bar{\beta }}_\omega \) is given in Sect. 6.4.5; the construction of the function \(\lambda _{\infty }\) is given in Sect. 6.4.7, while the construction of Z is given in Sect. 6.4.8; the proof of (5.1) is given in Sect. 6.6, and (5.2)–(5.3) can be proved along the same lines.

5.1 Proof of Theorem 1

The proof of the existence of the thermodynamic limit for the average of all the local functions of dimer configurations, Eq. (2.35), is postponed to Sect. 6.5, see in particular the comments after (6.149). The other statements in Theorem 1 follow easily from Proposition 2. In fact, if we rewrite Equation (5.1) by using Equations (4.17) and (4.19), we obtain Equation (2.36) with

$$\begin{aligned} {\bar{K}}_{\omega ,r}={{\hat{K}}}_{\omega ,r}\frac{1}{Z\sqrt{1-\tau ^2}}, \qquad {\bar{H}}_{\omega ,r}={{\hat{H}}}_{\omega ,r}\frac{\sqrt{{\bar{B}}}}{Z} \end{aligned}$$
(5.4)

which satisfy the desired properties, (2.40) and (2.44), thanks to the stated properties of \(Z, {{\hat{K}}}_{\omega ,r},{{\hat{H}}}_{\omega ,r}\) and of \({\bar{B}}\). Recall that \(\tau \) is defined in (4.18), and the critical exponent \(\nu \) is given in Equation (4.20): therefore, the fact that \(\nu =1+O(\lambda )\), see (2.39), follows from the definition there and from the fact that \(\lambda _\infty =O(\lambda )\). The other properties of \({\bar{p}}^\omega ,{\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \) stated in Theorem 1 are the same as those stated in items (1)–(2) of Proposition 2.

5.2 Proof of Theorem 2

The key ingredient in the proof of Theorem 2 is the analogue of (2.30)–(2.31) for the interacting case, namely formula (2.46). We start by discussing the proof of this formula which, as we shall see, is a direct consequence of the identities (5.2)–(5.3), and of the lattice Ward Identity (3.41). In fact, combining these three identities, we obtain:

(5.5)

where (with \(v_r\) as in (2.1))

$$\begin{aligned} {\mathcal {D}}_{\omega '}(p)=-i \sum _{r=2}^4{{\hat{K}}}_{\omega ',r}\,p\cdot v_r \end{aligned}$$

and, as before, \(0<{\mathfrak {c}}\le |p|,|k|,|k+p|\le 2{\mathfrak {c}}\). By comparing this equation with (4.24), and recalling that \({{\hat{v}}}(0)=1\) (see (4.12)) we get

$$\begin{aligned} Z(1-\tau )\sum _{\omega '=\pm }{\bar{D}}_{\omega '}(p){{\hat{G}}}^{(2,1)}_{R, \omega ',\omega }(k,p)=\sum _{\omega '=\pm }{\mathcal {D}}_{\omega '}(p) {{\hat{G}}}^{(2,1)}_{R,\omega ',\omega }(k,p) [1+O({\mathfrak {c}}^\theta )]\;. \end{aligned}$$
(5.6)

This implies that

$$\begin{aligned} -i\sum _{r=2}^4{{\hat{K}}}_{\omega ,r}(v_r)_1= Z(1-\tau ){\bar{\alpha }}_\omega ,\qquad -i\sum _{r=2}^4{{\hat{K}}}_{\omega ,r}(v_r)_2= Z(1-\tau ){\bar{\beta }}_\omega . \end{aligned}$$
(5.7)

In order to deduce (5.7) from (5.6), one can proceed as follows: by [26, Eq. C.24], we have that

$$\begin{aligned} G^{(2,1)}_{R,-\omega ,\omega }(k,p)=\tau {{\hat{v}}}(p) \frac{{\bar{D}}_\omega (p)}{{\bar{D}}_{-\omega }(p)} {{\hat{G}}}^{(2,1)}_{R,\omega ,\omega }(k,p). \end{aligned}$$
(5.8)

By plugging this identity into (5.6) we get (keeping the terms of dominant order as \(p\rightarrow 0\) only and using \({{\hat{v}}}(0)=1\)):

$$\begin{aligned} Z(1-\tau ^2){\bar{D}}_{\omega }(p){\bar{D}}_{-\omega }(p)={\mathcal {D}}_\omega (p) {\bar{D}}_{-\omega }(p)+ \tau {\bar{D}}_{\omega }(p){\mathcal {D}}_{-\omega }(p). \end{aligned}$$
(5.9)

Computing this formula at \(p_2=0, p_1\ne 0\) first, both for \(\omega =+\) and \(\omega =-\) and then repeating the computation for \(p_1=0, p_2\ne 0\), one gets a system of linear equations for the coefficients \(-i\sum _{r=2}^4{{\hat{K}}}_{\omega ,r}(v_r)_j\), with \(j=1,2, \omega =\pm \), whose solution is (5.7).

By replacing (5.4) into (5.7) and recalling that \(\nu =\frac{1-\tau }{1+\tau }\), cf. (4.20), we find

$$\begin{aligned} {\bar{K}}_{\omega ,2}+{\bar{K}}_{\omega ,3}=-i\sqrt{\nu }{\bar{\alpha }}_\omega ,\qquad {\bar{K}}_{\omega ,3}+{\bar{K}}_{\omega ,4}=-i\sqrt{\nu }{\bar{\beta }}_\omega . \end{aligned}$$
(5.10)

We claim that \(\sum _{r=1}^4{\bar{K}}_{\omega ,r}=0\) (we shall prove this fact in a moment): therefore, the first equation can be rewritten as \({\bar{K}}_{\omega ,1}+{\bar{K}}_{\omega ,4}=i\sqrt{\nu }{\bar{\alpha }}_\omega \). In terms of the ‘elementary steps’ s(xj) in direction \({{e}}_{j}\) centered at x, introduced before (2.30), the two equations in (5.10) become

$$\begin{aligned}&\sum _{e\in s(x,1)} \sigma _e {\bar{K}}_{\omega ,r(e)}=-i\sqrt{\nu }{\bar{\beta }}_\omega =-i\omega \sqrt{\nu }\Delta _1{\bar{\phi }}_\omega \end{aligned}$$
(5.11)
$$\begin{aligned}&\sum _{e\in s(x,2)} \sigma _e {\bar{K}}_{\omega ,r(e)}=i\sqrt{\nu }{\bar{\alpha }}_\omega = -i \omega \sqrt{\nu }\Delta _2{\bar{\phi }}_\omega , \end{aligned}$$
(5.12)

which are the two cases of Equation (2.46).

In order to complete the proof of (5.11)–(5.12), we need to prove that \(\sum _{r=1}^4{\bar{K}}_{\omega ,r}=0\), as claimed above. For this purpose, we consider (3.40), and combine it with (5.2)–(5.3), thus getting, if \(0<{\mathfrak {c}}\le |p|,|k|,|k+p|\le 2{\mathfrak {c}}\)

$$\begin{aligned} \sum _{r=1}^4\sum _{\omega '=\pm }{{\hat{K}}}_{\omega ',r} {{\hat{G}}}^{(2,1)}_{R,\omega ',\omega }(k,p)[1+O({\mathfrak {c}}^\theta )]{\mathop {=}\limits ^{ {\mathfrak {c}}\rightarrow 0}}{{\hat{G}}}^{(2)}_{R,\omega }(k+p)[1+O({\mathfrak {c}}^\theta )]. \end{aligned}$$
(5.13)

Using (5.8) to rewrite the left hand side and recalling that \({{\hat{v}}}(p)=1+O(p)\), this becomes

$$\begin{aligned}&{{\hat{G}}}^{(2,1)}_{R,\omega ,\omega }(k,p)\sum _{r=1}^4 \Big ({{\hat{K}}}_{\omega ,r}+\tau {{\hat{K}}}_{-\omega ,r} \frac{{\bar{D}}_\omega (p)}{{\bar{D}}_{-\omega }(p)}\Big )[1+O({\mathfrak {c}}^\theta )]=\nonumber \\&\quad ={{\hat{G}}}^{(2)}_{R,\omega }(k+p)[1+O({\mathfrak {c}}^\theta )]. \end{aligned}$$
(5.14)

Now, recalling that the magnitude of the correlation functions for small pk has the form given in (4.22) and (4.23), for this to hold in the limit \({\mathfrak {c}}\rightarrow 0\) we must have

$$\begin{aligned} \sum _{r=1}^4\Big ({{\hat{K}}}_{\omega ,r}+\tau {{\hat{K}}}_{-\omega ,r}\lim _{p_j\rightarrow 0}\frac{{\bar{D}}_\omega (p_j)}{{\bar{D}}_{-\omega }(p_j)}\Big )=0, \end{aligned}$$
(5.15)

for any sequence \(p_j\) along which the ratio \({\bar{D}}_\omega (p_j)/{\bar{D}}_{-\omega }(p_j)\) admits a limit. Note that, in general, the limit depends upon the chosen subsequence. For instance, if \(p_j=(s_j,0)\) with \(s_j\rightarrow 0\) then the limit is \(-\alpha _\omega ^2/|\alpha _\omega |^2\) while if \(p_j=(0,t_j)\) with \(t_j\rightarrow 0\) the limit is \(-\beta _\omega ^2/|\beta _\omega |^2\). On the other hand, these two values cannot be equal since we know that the ratio \(\alpha _\omega /\beta _\omega \) is not real (cf. (4.5)). Consequently, (5.15) implies that \(\sum _{r=1}^4{{\hat{K}}}_{\omega ,r}=0\) that, in light of (5.4), is equivalent to \(\sum _{r=1}^4 {\bar{K}}_{\omega ,r}=0\), as desired.

With the identities (5.11)–(5.12) at hand, we can easily prove (2.47), by repeating the analogue of the discussion leading, in the non-interacting case, to (2.33). We will be very sketchy since the analogous argument has been given in detail in [25] in the case of the model with weights \({\underline{t}}\equiv 1\). We start from the very definition of the covariance of the height difference:

$$\begin{aligned} {\mathbb {E}}_\lambda \left[ (h(\eta _1)-h(\eta _2));(h(\eta _3)-h(\eta _4)) \right] = \sum _{e\in C_{\eta _1\rightarrow \eta _2}} \sum _{e'\in C_{\eta _3\rightarrow \eta _4}}\sigma _e \sigma _{e'} {\mathbb {E}}_\lambda ({\mathbb {1}}_e;{\mathbb {1}}_{e'}), \end{aligned}$$
(5.16)

where \(C_{\eta _1\rightarrow \eta _2}\) and \(C_{\eta _3\rightarrow \eta _4}\) are two lattice paths connecting \(\eta _1\) with \(\eta _2\), and \(\eta _3\) with \(\eta _4\), respectively. For simplicity, we assume that \(\eta _1\) and \(\eta _2\) have the same parity, and similarly for \(\eta _3\) and \(\eta _4\): in this way, it is possible to choose the two paths \(C_{\eta _1\rightarrow \eta _2}\) and \(C_{\eta _3\rightarrow \eta _4}\) to be concatenations of ‘elementary steps’ s(xj) in directions \(\pm {{e}}_j\), see the discussion after (2.28) above. For simplicity, let us also assume that the mutual distances between the faces \(\eta _1,\dots \eta _4\) are all comparable, i.e.

$$\begin{aligned} 0< c<\frac{\min _{i\ne j}|\eta _i-\eta _j|}{\max _{i\ne j}|\eta _i-\eta _j|}. \end{aligned}$$
(5.17)

In this case, we choose the two paths \(C_{\eta _1\rightarrow \eta _2}\) and \(C_{\eta _3\rightarrow \eta _4}\) to be of length at most \(C \max _{i\ne j}|\eta _i-\eta _j|\) and to be at mutual distance \(C^{-1} \max _{i\ne j}|\eta _i-\eta _j|\), for some constant \(C=C(c)\).

We now insert (2.36) into (5.16) and, by repeating the discussion of [25, Section 3.2], we find that the dominant contribution comes from \({\bar{A}}_{r,r'}\) (the contribution from \({\bar{B}}_{r,r'}\) is sub-dominant due to the oscillating pre-factors):

$$\begin{aligned}&{\mathbb {E}}_\lambda \left[ (h(\eta _1)-h(\eta _2));(h(\eta _3)-h(\eta _4)) \right] = \end{aligned}$$
(5.18)
$$\begin{aligned}&\quad = \sum _{e\in C_{\eta _1\rightarrow \eta _2}} \sum _{e'\in C_{\eta _3\rightarrow \eta _4}}\sigma _e \sigma _{e'} {\bar{A}}_{r(e),r(e')}(x(e),x(e')) \end{aligned}$$
(5.19)
$$\begin{aligned}&\qquad +\, O\left( \frac{1}{\min _{i\ne j\le 4}|\eta _i-\eta _j|^{1/2}+1}\right) , \end{aligned}$$
(5.20)

where r(e) is the type of the edge e, x(e) is the coordinate of the black site of e. By using the explicit expression of \({\bar{A}}_{r,r'}\), (2.37), and by decomposing the two paths \(C_{\eta _1\rightarrow \eta _2}, C_{\eta _3\rightarrow \eta _4}\), into a sequence of elementary steps, we obtain (denoting the generic elementary step in \(C_{\eta _1\rightarrow \eta _2}\), resp. \(C_{\eta _3\rightarrow \eta _4}\), by s(xj), resp. \(s(x',j')\))

$$\begin{aligned} (5.18)= & {} \frac{1}{4\pi ^2}\sum _{\omega =\pm }\sum _{\begin{array}{c} s(x,j)\in C_{\eta _1\rightarrow \eta _2}\\ s(x',j')\in C_{\eta _3\rightarrow \eta _4} \end{array}} \sum _{\begin{array}{c} e\in s(x,j)\\ e'\in s(x',j') \end{array}} {\sigma _e\sigma _{e'}} \frac{{\bar{K}}_{\omega ,r(e)} {\bar{K}}_{\omega ,r(e')} }{({\bar{\phi }}_\omega (x-x'))^2}\nonumber \\&+\, O\left( \frac{1}{\min _{i\ne j\le 4}|\eta _i-\eta _j|^{1/2}+1}\right) . \end{aligned}$$
(5.21)

We now use (5.11)–(5.12) and the symmetry \({\bar{\phi }}_\omega ={\bar{\phi }}^*_{-\omega }\) to rewrite the dominant term in (5.21) as the Riemann sum approximating the following integral:

$$\begin{aligned} - \frac{\nu }{2\pi ^2}\mathfrak {R}\int _{{\bar{\phi }}_+(\eta _1)}^{{\bar{\phi }}_+(\eta _2)}dz\int _{{\bar{\phi }}_+(\eta _3)}^{{\bar{\phi }}_+(\eta _4)}dz' \frac{1 }{(z-z')^2} \end{aligned}$$
(5.22)

whose explicit evaluation gives the main term in the r.h.s. of (2.47). Putting together the error terms, we obtain the statement of Theorem 2, as desired.

In the case where (5.17) fails (e.g. when \(\eta _1=\eta _3\), \(\eta _2=\eta _4\) and (5.16) is just the variance of the height gradient), one chooses the paths \(C_{\eta _1\rightarrow \eta _2},C_{\eta _3\rightarrow \eta _4}\) to be “as well separated as possible” (cf. [25, Sec. 3.2]) and the rest of the argument works the same.

6 Renormalization Group Analysis

In this section we discuss the multiscale analysis of the dimer model and the comparison with the continuum model, which leads us to the results spelled out in Proposition 2.

We first make a few preliminary manipulations on the Grassmann integral, in order to set it up in a form appropriate for multiscale integration (Sect. 6.1). In the following sections, Sects. 6.2 to 6.6, we first give the iterative definition of the effective potentials, then explain how to bound the norm of their kernels, how to fix the dressed velocities and Fermi points in order to make the expansion convergent uniformly in the thermodynamic limit and, finally, how to adapt the bounds to the computation of the correlation functions. For a more detailed guidance on Sects. 6.2 to 6.6, see the end of Sect. 6.1.

From now on, \(C,C',\ldots ,\) and \(c,c',\ldots ,\) denote constants independent of L, whose specific values might change from line to line.

6.1 Preliminaries

As a preliminary step, we rewrite the quadratic part S of the action in (3.28) as a “dressed” term \(S_0\) plus a “counter-term” \(N=S-S_0\), whose role is to fix the location of the interacting Fermi points and Fermi velocities. Namely, letting as usual

$$\begin{aligned} {\hat{\psi }}^\pm _k=\sum _{x\in \Lambda }\psi ^\pm _xe^{\mp i k x}, \quad k\in {\mathcal {P}}({\varvec{\theta }}), \qquad \psi ^\pm _x=\frac{1}{L^2}\sum _{k\in {\mathcal {P}}({\varvec{\theta }})}e^{\pm ikx}{\hat{\psi }}^\pm _k, \end{aligned}$$
(6.1)

we write:

$$\begin{aligned} S(\psi )=-L^{-2}\sum _{k\in {\mathcal {P}}({\varvec{\theta }})} \mu (k){\hat{\psi }}^+_k{\hat{\psi }}^-_k\equiv S_0(\psi )+N(\psi ), \end{aligned}$$
(6.2)

where \(S_0(\psi )=-L^{-2}\sum _{k\in {\mathcal {P}}({\varvec{\theta }})} \mu _0 (k){\hat{\psi }}^+_k{\hat{\psi }}^-_k\), with

$$\begin{aligned} \mu _0(k)=\mu (k)+\sum _{\omega =\pm }{\bar{\chi }}_0(k-{\bar{p}}^\omega ) \left[ -\mu ({\bar{p}}^\omega )+a_\omega (k_1-{\bar{p}}^\omega _1) +b_\omega (k_2-{\bar{p}}^\omega _2)\right] . \end{aligned}$$
(6.3)

In this equation:

  1. (1)

    \({\bar{p}}^\omega ={\bar{p}}^\omega (\lambda )\), with \(\omega =\pm \), are points in \([-\pi ,\pi ]^2\), such that \({\bar{p}}^++{\bar{p}}^-=(\pi ,\pi )\), and they will be fixed via the multiscale construction. A posteriori they can be interpreted as “dressed Fermi points”; they are the same functions appearing in Theorem 1.

  2. (2)

    \(a_\omega =a_\omega (\lambda )\in {\mathbb {C}}\) and \(b_\omega =b_\omega (\lambda )\in {\mathbb {C}}\) are such that \(a_\omega =-a^*_{-\omega }\) and \(b_\omega =-b^*_{-\omega }\); they will also be fixed via the multiscale construction. A posteriori, their choice fixes the “dressed Fermi velocities” via the following relations:

    $$\begin{aligned} \partial _{p_1}\mu _0({\bar{p}}^\omega )= & {} \partial _{p_1}\mu ({\bar{p}}^\omega )+a_\omega =:{\bar{\alpha }}_\omega , \end{aligned}$$
    (6.4)
    $$\begin{aligned} \partial _{p_2}\mu _0({\bar{p}}^\omega )= & {} \partial _{p_2}\mu ({\bar{p}}^\omega ) +b_\omega =:{\bar{\beta }}_\omega , \end{aligned}$$
    (6.5)

    where \({\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \) are the same functions appearing in Theorem 1.

  3. (3)

    the function \({\bar{\chi }}_0\) is defined as: \({\bar{\chi }}_0(k')={\bar{\chi }}(|{\mathcal {M}}^{-1}k'|)\), where: (1) \({\mathcal {M}}\) is the same matrix as (4.1), with \({\bar{\alpha }}^{1}\) and \({\bar{\alpha }}^2\) (resp. \({\bar{\beta }}^{1}\) and \({\bar{\beta }}^2\)) the real and imaginary parts of \({\bar{\alpha }}_+\) (resp. \({\bar{\beta }}_+\)); (2) \({\bar{\chi }}:{\mathbb {R}}^+\rightarrow [0,1]\) is a \(C^\infty \) cut-off function in the Gevrey class of order 2 (see [25, Appendix C]) that is equal to 1, if its argument is smaller than \(c_0/2\), and equal to 0, if its argument is larger than \(c_0\); here \(c_0\) is a small enough constant, such that in particular the support of \({\bar{\chi }}_0(\cdot -{\bar{p}}^+)\) is disjoint from the support of \({\bar{\chi }}_0(\cdot -{\bar{p}}^-)\). For later reference, we also let for h a negative integer

    $$\begin{aligned} {\bar{\chi }}_h(k'):={\bar{\chi }}_0(2^{-h}k') . \end{aligned}$$
    (6.6)

From the properties just stated of \({\bar{p}}^\omega ,a_\omega ,b_\omega \) and \({\bar{\chi }}(\cdot )\), we see that

$$\begin{aligned} \mu _0((\pi ,\pi )-k)=\mu _0^*(k). \end{aligned}$$
(6.7)

In the integration over \(\psi \) in (3.28), the Fourier modes k that are the closest to the zeros of \(\mu _0(\cdot )\) play a somewhat special role, so they have to be treated separately, at the very last step of the multi-scale procedure (see Sect. 6.5). Namely, given \({\varvec{\theta }}\in \{0,1\}^2\), let \(k^\pm _{{\varvec{\theta }}}\) be the values of \(k\in {\mathcal {P}}({\varvec{\theta }})\) that are closest to \({\bar{p}}^\pm \) and note that \(k^+_{\varvec{\theta }}=(\pi ,\pi )-k^-_{\varvec{\theta }}\) [If there is more than one momentum at minimal distance from \({\bar{p}}^\pm \) (there are at most four), any arbitrary choice will work]. Next, we decompose the quadratic action \(S_0(\psi )\) as a sum of a term depending only on \(k^\pm _{\varvec{\theta }}\) plus a term depending only on the modes in

$$\begin{aligned} {\mathcal {P}}'({\varvec{\theta }}):={\mathcal {P}}({\varvec{\theta }}){\setminus } \{k^+_{\varvec{\theta }},k^-_{\varvec{\theta }}\}, \end{aligned}$$

and we rewrite (3.28) as

$$\begin{aligned} e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,\phi )}= & {} \int D\psi \,e^{-L^{-2}\sum _{\omega =\pm }\mu _0(k_{\varvec{\theta }}^\omega ){\hat{\psi }}^+_{k_{\varvec{\theta }}^\omega } {\hat{\psi }}^-_{k_{\varvec{\theta }}^\omega }}\nonumber \\&\times e^{-L^{-2}\sum _{k\in {\mathcal {P}}'({\varvec{\theta }})} \mu _0(k){\hat{\psi }}^+_{k}\psi ^-_{k}+N(\psi )+V(\psi ,A)+(\psi ,\phi )}. \end{aligned}$$
(6.8)

We multiply and divide by

$$\begin{aligned} e^{L^2 E^{(0)}}:=\prod _{k\in {\mathcal {P}}'({\varvec{\theta }})} \mu _0(k), \end{aligned}$$
(6.9)

(the product is non-zero since we singled out the possibly zero modes \(k^\pm _{{\varvec{\theta }}}\)) and, letting

$$\begin{aligned} {\hat{\Psi }}_\omega ^\pm :={\hat{\psi }}^\pm _{k_{\varvec{\theta }}^\omega }, \end{aligned}$$
(6.10)

we rewrite the generating function as

$$\begin{aligned} e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,\phi )}=\int D{\hat{\Psi }} e^{-L^{-2}\sum _{\omega =\pm }\mu _0(k_{\varvec{\theta }}^\omega ){\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_ {\omega } +{\mathbb {W}}_L^{({\varvec{\theta }})}(A,\phi ,\Psi )}. \end{aligned}$$
(6.11)

Here

$$\begin{aligned} \Psi ^\pm _x=\frac{1}{L^2}\sum _{\omega =\pm }e^{\pm i k_{\varvec{\theta }}^\omega x}{\hat{\Psi }}^\pm _{\omega }, \quad \int D{\hat{\Psi }}\prod _{\omega =\pm } {\hat{\Psi }}^-_{\omega }{\hat{\Psi }}^+_{\omega }= L^4, \end{aligned}$$
(6.12)

(the \(L^4\) factor comes from the fact that (3.10) translates in Fourier space into \( \int D\psi \prod _{k\in {\mathcal {P}}}{({\varvec{\theta }})}[L^{-2}{{\hat{\psi }}^-_k{\hat{\psi }}^+_k}]=1\)) and

$$\begin{aligned} e^{{\mathbb {W}}_L^{({\varvec{\theta }})}(A,\phi ,\Psi )}:=e^{L^2 E^{(0)}} \int P_{g_{0}}(d\psi ) e^{N(\Psi +\psi )+V(\Psi +\psi ,A)+(\Psi +\psi ,\phi )}, \end{aligned}$$
(6.13)

with \(P_{g_{0}}\) the Grassmann Gaussian integration (see footnote 6) with propagator

$$\begin{aligned} g_{0}(x,y)=L^{-2}\sum _{k\in {\mathcal {P}}'({\varvec{\theta }})} \frac{e^{-ik(x-y)}}{\mu _0(k)}. \end{aligned}$$
(6.14)

From this point, we proceed as follows. First, we perform in a multi-scale way the integration over the Grassmann variables \(\psi \), i.e. over the Fourier modes except \(k^\pm _{{\varvec{\theta }}}\): the inductive integration procedure, including the definition of the running coupling constants (RCC), is described in Sect. 6.2; the outcome of the construction can be conveniently expressed in terms of a Gallavotti-Nicolò tree expansion, similar to the one described in [25, Section 6.2]. The main definitions (and the main differences compared to the case treated in [25]) are summarized in Sect. 6.3; in the same section, we also state the bounds satisfied by the kernels of the effective potential, see Proposition 3, under the assumption that the RCC are uniformly bounded in the infrared, see condition (6.66). The proof that the RCC remain in fact bounded under the iterations of the renormalization group map is given in Sect. 6.4; the flow of the RCC can be controlled only if their initial data are properly fixed: as shown there, the choice of the initial data fixes the dressed Fermi points \({\bar{p}}^\omega \) and the dressed Fermi velocities \({\bar{\alpha }}_\omega \), \({\bar{\beta }}_\omega \), as anticipated after (6.3). In Sect. 6.5 we describe the integration of the last two modes and prove the existence of the thermodynamic limit for the correlation functions, with explicit bounds on the speed of convergence as \(L\rightarrow \infty \). Finally, in Sect. 6.6, we compute the fine asymptotics of the correlations functions, via a comparison of the tree expansion of the dimer model with that of the continuum model of Sect. 4, and complete the proof of Proposition 2.

6.2 Multi-scale analysis

In this section we describe the multi-scale computation of \({\mathbb {W}}_L^{({\varvec{\theta }})}(A,\phi , \Psi )\) defined in (6.13). We consider explicitly only the case \(\phi =0\); the general case can be treated analogously but we will not belabor the details in this paper.

The procedure is based on a systematic use of the ‘addition principle’ for Gaussian Grassmann integrals, namely the following property [21, Sec. 4]: if \(P_{g}(d\psi )\) is the Grassmann Gaussian integration with propagator g and \(g=g_1+g_2\) then

$$\begin{aligned} \int P_{g}(d\psi ) F(\psi )=\int P_{g_1}(d\psi _1) P_{g_2}(d\psi _2) F(\psi _1+\psi _2). \end{aligned}$$
(6.15)

We apply this formula to \(P_{g_{0}}\), in connection with the following decomposition of the propagator \(g_{0}(x,y)\):

$$\begin{aligned} g_{0}(x,y)=g^{(0)}(x,y)+\sum _{\omega =\pm } e^{-i{\bar{p}}^\omega (x-y)}g_{\omega }^{(\le -1)}(x,y) \end{aligned}$$
(6.16)

where

$$\begin{aligned} g^{(0)}(x,y)= L^{-2}\sum _{k\in {\mathcal {P}}'({\varvec{\theta }})}e^{-i k (x-y)}\frac{1-{\bar{\chi }}_{-1}(k-{\bar{p}}^+)-{\bar{\chi }}_{-1}(k-{\bar{p}}^-)}{\mu _0(k)} \end{aligned}$$
(6.17)

and, if \({\mathcal {P}}_\omega '({\varvec{\theta }})=\{k': k'+{\bar{p}}^\omega \in {\mathcal {P}}'({\varvec{\theta }})\}\),

$$\begin{aligned} g^{(\le -1)}_{\omega }(x,y)=L^{-2}\sum _{k'\in {\mathcal {P}}_\omega '({\varvec{\theta }})}e^{-ik'(x-y)}\frac{{\bar{\chi }}_{-1}(k')}{\mu _0(k'+{\bar{p}}^\omega )}. \end{aligned}$$
(6.18)

By using the decomposition (6.16) and (6.15), we rewrite (6.13) as

$$\begin{aligned}&e^{{\mathbb {W}}_L^{({\varvec{\theta }})}(A,0,\Psi )}= e^{L^2 E^{(0)}}\int P_{{(\le -1)}}(d\psi ^{(\le -1)})\times \nonumber \\&\quad \times \int P_{(0)}(d\psi ^{(0)}) e^{N(\psi ^{(0)}+\psi ^{(\le -1)}+\Psi )+V(\psi ^{(0)}+\psi ^{(\le -1)}+\Psi ,A)}, \end{aligned}$$
(6.19)

where \(\psi ^{(0)}+\psi ^{(\le -1)}+\Psi \) is a shorthand notation for

$$\begin{aligned} \{\psi ^{(0)\pm }_x+\sum _{\omega }e^{\pm i{\bar{p}}^\omega x}\varphi ^\pm _{x,\omega }\}_{x\in \Lambda }, \quad \varphi ^\pm _{x,\omega }:=\psi ^{(\le -1)\pm }_{x,\omega }+L^{-2}e^{\pm i (k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega )x}{\hat{\Psi }}^\pm _\omega .\quad \end{aligned}$$
(6.20)

\(P_{(0)}\) is the Grassmann Gaussian measure with propagator \(g^{(0)}(x,y)\), while \(P_{(\le -1)}\) is the Grassmann Gaussian measure with propagator

$$\begin{aligned} \delta _{\omega ,\omega '}g^{(\le -1)}_{\omega }(x,y)=\int P_{(\le -1)}(d\psi ) \psi ^{(\le -1)-}_{x,\omega }\psi ^{(\le -1)+}_{y,\omega '}. \end{aligned}$$

Since the cutoff function \({\bar{\chi }}_{-1}\) in (6.18) is a Gevrey function of order 2, the propagator \(g^{(0)}\) has stretched-exponential decay at large distances:

$$\begin{aligned} |g^{(0)}(x,y)|\le C e^{-\kappa \sqrt{|x-y|}}, \end{aligned}$$
(6.21)

for suitable L-independent constants \(C,\kappa >0\), if \(|x-y|\) is the distance on the torus \(\Lambda \). This is seen by writing \(g^{(0)}\) via the Poisson summation formula as a sum of Fourier integrals, as in [25, App. A]; each integral decays in the desired way because it is the Fourier transform of a Gevrey function [39].

Next, we denote by \(V^{(0)}(\cdot ,J)\) the combination \(N(\cdot )+V(\cdot ,A)\), re-expressed in terms of the variables \(J=\{J_{x,r}\}_{x\in \Lambda ,\, 1\le r\le 4}\), instead of A: here, if b is the bond of type r and black site x, we let \(J_{x,r}:=e^{A_b}-1\). The result of the integration over \(\psi ^{(0)}\) is rewritten in exponential form:

$$\begin{aligned} e^{L^2 E^{(0)}}\int P_{(0)}(d\psi ^{(0)}) e^{V^{(0)}(\psi ^{(0)}+\varphi ,J)}= e^{L^2 E^{(-1)} + S^{(-1)}(J)+V^{(-1)}(\varphi ,J)}, \end{aligned}$$
(6.22)

where [21, Sec. 4]

$$\begin{aligned}&L^2 (E^{(-1)}-E^{(0)}) +S^{(-1)}(J)+ V^{(-1)}(\varphi ,J)=\nonumber \\&\qquad =\sum _{n\ge 1}\frac{1}{n!}{\mathcal {E}}^T_0(\underbrace{V^{(0)}(\psi ^{(0)}+\varphi ,J);\cdots ; V^{(0)}(\psi ^{(0)}+\varphi ,J))}_{n\ \mathrm{times}}, \end{aligned}$$
(6.23)

with \({\mathcal {E}}^T_0\) the truncated expectationFootnote 10 w.r.t. the Grassmann Gaussian integration \(P_{(0)}(d\psi ^{(0)})\), and \(E^{(-1)}, S^{(-1)}(\cdot )\) are fixed by the condition \(S^{(-1)}(0)=0\), \(V^{(-1)}(0,J)=0\). The series in the r.h.s. is absolutely summable, for \(\lambda \) sufficiently small (independently of L), see [21, Sec. 4.2]. The reason is that the propagator \(g^{(0)}\) has a fast decay in space, uniformly in L, as in (6.21).

The effective potential on scale \(-1\) can be represented as in the following formula (which is a definition of the kernels \(W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}\)):

$$\begin{aligned}&V^{(-1)}(\varphi ,J)= \sum _{\begin{array}{c} n,m\ge 0:\\ n\ \mathrm{even}, \ n\ge 2 \end{array}}\sum _{{\underline{x}},\, {\underline{y}},\, {\underline{\omega }},\,{\underline{r}}}W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}({\underline{x}}, {\underline{y}})\nonumber \\&\qquad \times \, \varphi ^+_{x_1,\omega _1}\varphi ^-_{x_2,\omega _2}\cdots \varphi ^+_{x_{n-1},\omega _{n-1}}\varphi ^-_{x_n,\omega _n} J_{y_1,r_1}\cdots J_{y_m,r_m}, \end{aligned}$$
(6.24)

where: \({\underline{x}}=(x_1,\ldots , x_n)\in \Lambda ^n,{\underline{y}}\in \Lambda ^m,{\underline{\omega }}\in \{-1,+1\}^n,{\underline{r}}\in \{1,\dots ,4\}^m\); the Grassmann variables \(\varphi ^\pm _{x,\omega }\) were defined in (6.20). Moreover, the kernels can be written as

$$\begin{aligned} W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}({\underline{x}}, {\underline{y}})={{\tilde{W}}}^{(-1)}_{n,m;{\underline{r}}}({\underline{x}}, {\underline{y}})\,\exp \{i\sum _{j=1}^n(-1)^{j-1}{\bar{p}}^{\omega _j}x_j\}, \end{aligned}$$
(6.25)

with \({{\tilde{W}}}^{(-1)}_{n,m;{\underline{r}}}({\underline{x}}, {\underline{y}})\) a function that isFootnote 11 independent of \({\underline{\omega }}\), translationally invariant, periodic of period L in \(y_i\), and \({\varvec{\theta }}\)-periodic of period L in \(x_i\) (here we say that, e.g., a function is (0, 1)-periodic if it is periodic in the first coordinate and anti-periodic in the second, and similarly for the other cases). The kernels \(W^{(-1)}\) are not uniquely defined by (6.24); due to the anti-commutation of Grassmann variables and to the fact that \(J_{y,r}\) are ordinary commuting variables, one can assume that \(W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}({\underline{x}}, {\underline{y}})\) are symmetric under permutations of the indices \((y_1,r_1), \ldots , (y_m,r_m)\), and anti-symmetric under permutations of the indices \(\{(x_{2i},\omega _{2i})\}_{1\le i\le n/2}\) and of the indices \(\{(x_{2i-1},\omega _{2i-1})\}_{1\le i\le n/2}\). An analogous representation is valid for \(S^{(-1)}(\cdot )\), and we denote its kernels by \(W^{(-1)}_{0,m;{\underline{\omega }},{\underline{r}}}({\underline{y}})\).

There is an equivalent expression for \(V^{(-1)}\) in Fourier space. We use the following convention for the Fourier transforms of the fields \(\psi ,J\):

$$\begin{aligned} \varphi ^\pm _{x,\omega }=L^{-2}\sum _{k\in {\mathcal {P}}_\omega ({\varvec{\theta }})}e^{\pm i k\cdot x}{\hat{\varphi }}^\pm _{k,\omega },\qquad J_{x,r}=L^{-2}\sum _{p\in {\mathcal {P}}(\mathbf{0})} {{\hat{J}}}_{p,r} e^{-ipx}, \end{aligned}$$

where \({\mathcal {P}}_\omega ({\varvec{\theta }}):=\{k: k+{\bar{p}}^\omega \in {\mathcal {P}}({\varvec{\theta }})\}\). The reason why \(k\in {\mathcal {P}}_\omega ({\varvec{\theta }})\) (and not in \({\mathcal {P}}({\varvec{\theta }})\) as in (6.1)) is that the combination \(e^{\pm i {\bar{p}}^\omega x}\varphi ^\pm _{x,\omega }\) is \({\varvec{\theta }}\)-periodic, and not \(\varphi ^\pm _{x,\omega }\) itself. This sum includes also the momenta \(k=k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega \), \(\omega =\pm \). Of course, recalling that the only non-zero modes of \(\psi _{x,\omega }^\pm \) (resp. \(\Psi ^\pm _{x,\omega }\)) are in \({\mathcal {P}}'_\omega ({\varvec{\theta }})\) (resp. \(\{k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega \}_{\omega =\pm }\)), we have that

$$\begin{aligned} {\hat{\varphi }}_{k,\omega }^\pm ={\left\{ \begin{array}{ll} {\hat{\psi }}^\pm _{k,\omega }, \quad {\mathrm{if}}\ k\in {\mathcal {P}}'_\omega ({\varvec{\theta }}),\\ {\hat{\Psi }}^\pm _{\omega }, \quad {\mathrm{if}}\ k\not \in {\mathcal {P}}'_\omega ({\varvec{\theta }}). \end{array}\right. } \end{aligned}$$

Then, (6.24) becomes

$$\begin{aligned}&V^{(-1)}(\varphi ,J)= \sum _{\begin{array}{c} n,m\ge 0:\\ n\ {\mathrm{even}}, \ n\ge 2 \end{array}}\quad L^{-2(n+m)}\quad \sum _{{\underline{k}},\, {\underline{p}},\, {\underline{\omega }},\,{\underline{r}}} \quad {{\hat{W}}}^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}(k_2,\ldots ,k_n,p_1,\ldots ,p_m) \nonumber \\&\qquad \times \, {\hat{\varphi }}^+_{k_1,\omega _1}{\hat{\varphi }}^-_{k_2,\omega _2} \cdots {\hat{\varphi }}^+_{k_{n-1},\omega _{n-1}}{\hat{\varphi }}^-_{k_n,\omega _n} {{\hat{J}}}_{p_1,r_1}\cdots {{\hat{J}}}_{p_m,r_m}\delta _{{\underline{\omega }}}({\underline{k}}, {\underline{p}}), \end{aligned}$$
(6.26)

where \({\underline{k}}=(k_1,\ldots , k_n)\), with \(k_i\in {\mathcal {P}}_{\omega _i}({\varvec{\theta }})\), \({\underline{p}}=(p_1,\ldots p_m)\in \big [{\mathcal {P}}(\mathbf{0})\big ]^m\) and

(6.27)

is the periodized Kronecker delta enforcing momentum conservation. Also, \({{\hat{W}}}^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}(k_2,\ldots ,k_n,p_1,\ldots ,p_m)\) is just the Fourier transform of \({{\tilde{W}}}^{(-1)}_{n,m;{\underline{r}}}\), computed at momenta \(k_2+{\bar{p}}^{\omega _2},\ldots , k_n+{\bar{p}}^{\omega _n},p_1,\ldots ,p_m\) (it depends only on \(n+m-1\) momenta, due to translation invariance of \({{\tilde{W}}}^{(-1)}_{n,m;{\underline{r}}}\) in real space).

Using the Battle-Brydges-Federbush-Kennedy (BBFK) determinant formula and the Gram-Hadamard bound [21, Sec. 4.2] for the truncated expectation in (6.23), we find that \(E^{(-1)}\), and \(W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}({\underline{x}}, {\underline{y}})\) are absolutely convergent series and real analytic functions of

$$\begin{aligned} (\nu _{0,\omega },a_{0,\omega },b_{0,\omega },\lambda _0), \end{aligned}$$
(6.28)

for \(\max \{|\nu _{0,\omega }|,|a_{0,\omega }|,|b_{0,\omega }|,|\lambda _0|\}\le \varepsilon \) with \(\varepsilon \) sufficiently small, where we denoted (for uniformity of notation with the running coupling constants \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega },\lambda _h\), to be introduced below):

$$\begin{aligned} \nu _{0,\omega }:=-\mu ({\bar{p}}^\omega ),\quad a_{0,\omega }:=a_\omega ,\quad b_{0,\omega }:=b_\omega ,\quad \lambda _0:=\lambda . \end{aligned}$$
(6.29)

Moreover, \(|E^{(-1)}|\le C\varepsilon \) and, using also the exponential decay of the bare potential, (3.19), we find that

$$\begin{aligned} \Vert W^{(-1)}_{n,m}\Vert _{\kappa ,-1}\le C^{n+m}\varepsilon ^{\max \{1,cn\}}, \end{aligned}$$
(6.30)

for suitable constants \(\kappa ,C,c>0\) independent of the system size. Here

$$\begin{aligned} \Vert W^{(-1)}_{n,m}\Vert _{\kappa ,-1}:=L^{-2}\sup _{{\underline{\omega }}, {\underline{r}}}\sum _{{\underline{x}}, {\underline{y}}}|W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}({\underline{x}}, {\underline{y}})|e^{\kappa \sqrt{2^{-1}d({\underline{x}},{\underline{y}})}}, \end{aligned}$$
(6.31)

and \(d(x_1,\ldots ,x_l)\) is the length of the shortest tree on the torus \(\Lambda \) connecting the l points in \((x_1,\ldots ,x_l)\). The choice of the stretched-exponential weight in (6.31) is related to the stretched-exponential decay of the propagator, see (6.21). For technical details about the proof (6.31), or, better, of its analogue in a similar context, the reader can consult, e.g., [22, Section III.A and Eq. (3.19)].

Remark 8

The fact that the kernels \(W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}\) are absolutely convergent series of \((\nu _{0,\omega },\)\(a_{0,\omega },b_{0,\omega },\lambda _0)\), that each term in the expansion admits a limit as \(L\rightarrow \infty \) (as one can check by inspection) and that they satisfy uniform bounds as \(L\rightarrow \infty \), see (6.30), implies that their infinite volume limits exist and satisfy the same bounds. For later reference, the infinite volume limit of \(W^{(-1)}_{n,m;{\underline{\omega }},{\underline{r}}}\) will be denoted by \(W^{(-1),\infty }_{n,m;{\underline{\omega }},{\underline{r}}}\), and similarly for its Fourier transform.

After this first integration step, we still need to integrate \(\psi ^{(\le -1)}\) out, see (6.19). Let us first informally explain how this is done, before giving the precise inductive procedure in Sects. 6.2.16.2.3. The idea is to repeat the same procedure as above: we rewrite (via the addition principle) \(\psi ^{(\le -1)}_{\omega } = \psi ^{(-1)}_\omega +\psi ^{(\le -2)}_{\omega }\), where \(\psi ^{(-1)}_\omega \) (resp. \(\psi ^{(\le -2)}_\omega \)) is a Grassmann field with propagator supported, in momentum space, on momenta \(k'\in {\mathcal {P}}'_\omega ({\varvec{\theta }})\) with \(|k'|\sim 2^{-1}\) (resp. \(|k'|\lesssim 2^{-2}\)); we integrate \(\psi ^{(-1)}_\omega \) out; we exponentiate the result of the integration, thus defining the effective potential on scale \(-2\), in analogy with (6.22)–(6.23); and so on. One after the other, we integrate the fields \(\psi ^{(-2)}, \ldots , \psi ^{(h+1)}\) out, define the effective potential \(V^{(h)}\) on scale h (which involves fields \(\psi ^{(\le h)}\) with momenta \(k'\in {\mathcal {P}}'_\omega ({\varvec{\theta }})\) that belong to the support of \({\bar{\chi }}_{h}(\cdot )\) (cf. (6.6)), and continue until we reach the ‘last scale’, \(h_L\), fixed by the finite volume L, which induces a natural infrared cut-off. More precisely, \(h_L\) is fixed as the smallest (in absolute value) negative integer h such that the support of \({\bar{\chi }}_{h}(\cdot )\) has empty intersection with \({\mathcal {P}}'_\omega ({\varvec{\theta }})\). Note that, since all momenta in \({\mathcal {P}}'_\omega ({\varvec{\theta }})\) are at distance at least \(\pi /L\) from \({\bar{p}}^\omega \), we have \(h_L\sim -\log _2L\) for L large. The result of the integration of the Grassmann fields \(\psi ^{(\le h_L)}\) gives the generating function \({\mathbb {W}}_L^{({\varvec{\theta }})}(A,0,\Psi )\), as desired.

In order for the bounds on the generating function to be uniform in L, we need to improve the procedure roughly described here: at each step, before integrating the field on the next scale, we actually need to isolate and re-sum a certain selection of potentially dangerous contributions to the effective potential, the so-called marginal and relevant terms. We refer, e.g., to [25, Sec. 5], see in particular [25, Section 5.2.2] for a dimensional classification of the divergent terms arising in a ‘naive’ multiscale scheme. As discussed there, see [25, Eq. (5.8)] and following lines, the scaling dimension of the kernels with n external fields of type \(\psi \) and m external fields of type J is \(2-n/2-m\); in the renormalization group jargon, positive scaling dimension (that is, \(2-n/2-m>0\)\(\Leftrightarrow \)\((n,m)=(2,0)\)) corresponds to relevant contributions, vanishing scaling dimension (that is, \(2-n/2-m=0\)\(\Leftrightarrow \)\((n,m)=(4,0),(2,1)\)) corresponds to marginal contributions, and negative scaling dimension corresponds to irrelevant ones. In order to cure the potential divergences associated with the terms with \((n,m)=(2,0),(4,0),(2,1)\), at each step of the multiscale construction we properly ‘localize’ and re-sum these terms, via an iterative procedure that we now describe.

6.2.1 The inductive statement

Let us inductively assume that the fields \(\psi ^{(0)}, \psi ^{(-1)},\)\(\ldots ,\psi ^{(h+1)}\), \(h\ge h_L\), have been integrated out, and that after their integration the generating function has the following structure, analogous to the one at scales \(0,-1\):

$$\begin{aligned}&e^{-L^{-2}\sum _\omega \mu _0(k^\omega _{\varvec{\theta }}){\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_{\omega } +{\mathbb {W}}_L^{({\varvec{\theta }})}(A,0,\Psi )}=e^{L^2 E^{(h)}+S^{(h)}(J)} \nonumber \\&\qquad \times e^{-L^{-2}Z_h\sum _\omega \mu _{h,\omega }(k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega ){\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_{\omega }}\int P_{(\le h)}(d\psi )e^{V^{(h)}(\sqrt{Z_h}(\psi +\Psi ),J)}, \end{aligned}$$
(6.32)

for suitable real constants \(E^{(h)}\), \(Z_h\), and suitable ‘effective potentials’ \(S^{(h)}(J)\), \(V^{(h)}(\varphi ,J)\), to be defined inductively below, and fixed in such a way that \(V^{(h)}(0,J)=S^{(h)}(0)=0\). In the second line,

$$\begin{aligned} \mu _{h,\omega }(k):={\bar{D}}_\omega (k)+r_\omega (k)/{Z_h}, \end{aligned}$$

where

$$\begin{aligned} {\bar{D}}_\omega (k)={\bar{\alpha }}_\omega k_1+{\bar{\beta }}_\omega k_2 \end{aligned}$$

and

$$\begin{aligned} r_\omega (k)=\mu (k+{\bar{p}}^\omega ) -\mu ({\bar{p}}^\omega )-\partial _{k_1}\mu ({\bar{p}}^\omega )k_1 -\partial _{k_2}\mu ({\bar{p}}^\omega )k_2 \end{aligned}$$
(6.33)

is a remainder of order \(O(k^2)\) for k small. Finally, \(P_{(\le h)}(d\psi )\) is the Grassmann Gaussian integration with propagator (diagonal in the index \(\omega \))

$$\begin{aligned} \frac{1}{Z_h}g^{(\le h)}_{\omega }(x,y)=\frac{1}{Z_h}\frac{1}{L^2} \sum _{k\in {\mathcal {P}}_\omega '({\varvec{\theta }})}e^{-ik(x-y)} \frac{{\bar{\chi }}_{h}(k)}{\mu _{h,\omega }(k)}. \end{aligned}$$
(6.34)

We will also prove inductively that:

  1. (1)

    \(V^{(h)}(\varphi ,J)\) has the same structure as (6.26), with the upper index \((-1)\) in the kernels replaced by (h);

  2. (2)

    the kernels of \(V^{(h)}(\varphi ,J)\) satisfy the following symmetry:

    $$\begin{aligned} {\hat{W}}^{(h)}_{n,m;-{\underline{\omega }},{\underline{r}}}({\underline{k}},{\underline{p}})= \Bigl [{{\hat{W}}}^{(h)}_{n,m;{\underline{\omega }},{\underline{r}}}(-{\underline{k}},-{\underline{p}})\Bigr ]^*. \end{aligned}$$
    (6.35)

Remark 9

It is important to emphasize right away that we will view the kernels \(W^{(h)}_{n,m;{\underline{\omega }},{\underline{r}}}\), \(h\le -2\), as functions of:

  1. (i)

    a sequence of running coupling constants

    $$\begin{aligned} \{\lambda _{h'},\nu _{h',\omega },a_{h',\omega },b_{h',\omega },Y_{h',r,(\omega ,\omega ')}\}_{h< h'\le -1}. \end{aligned}$$
  2. (ii)

    a sequence of single-scale propagators \(\{g^{(h')}_\omega /Z_{h'-1}\}_{h< h'\le -1}\), of the form

    $$\begin{aligned} \frac{1}{Z_{h-1}}g^{(h)}_\omega (x,y):=\frac{1}{L^2}\sum _{k\in {\mathcal {P}}_\omega '({\varvec{\theta }})}e^{-ik(x-y)}\frac{f_{h}(k)}{{{\tilde{Z}}}_{h-1}(k)\bar{D}_\omega (k)+r_\omega (k)}, \end{aligned}$$
    (6.36)

    where \(f_{h}(k)={\bar{\chi }}_h(k)-{\bar{\chi }}_{h-1}(k)\) and

    $$\begin{aligned} {{\tilde{Z}}}_{h-1}(k)=Z_{h-1}{\bar{\chi }}_h(k)+Z_h(1-{\bar{\chi }}_{h}(k)); \end{aligned}$$
  3. (iii)

    the irrelevant part of \(V^{(-1)}\), denoted by \({\mathcal {R}} V^{(-1)}\).

The actual values of the running coupling constants (RCCs) will be defined via an inductive procedure in Sects. 6.2.26.2.3, the outcome of which is the beta function equation (6.57). In Sect. 6.4, we will show that there is only one specific choice of the initial data \((\nu _{0,\omega },a_{0,\omega },b_{0,\omega })\), which we will determine via a fixed point argument, guaranteeing that the flow of RCCs is uniformly bounded for all \(h\le 0\). For the fixed point argument itself, it is convenient to allow the beta function, as well as the kernels of the effective potential, to be computed at values of the RCCs different from the final, ‘correct’, ones. This is what we mean by saying that \(W^{(h)}_{n,m;{\underline{\omega }},{\underline{r}}}\) will be thought of as functions of the RCCs: we will allow ourselves to think of the RCCs as independent variables, which can be varied freely, as long as they remain sufficiently small; similarly for the dependences on \(g^{(h)}_\omega /Z_{h-1}\) and \({\mathcal {R}} V^{(-1)}\) mentioned in items (ii)–(iii): for certain manipulations discussed below, we will allow ourselves to modify the definition of the kernels by modifying the form of the single-scale propagators or of the kernel of the irrelevant part at scale \(-1\), keeping the rest of the iterative definition unchanged.

6.2.2 The inductive statement for \(h=-1\)

The representation (6.32) with (6.34)–(6.33) is valid at the initial step, \(h=-1\), with \(Z_{-1}=1\). To see this, one needs to use that, if k belongs to the support of \({\bar{\chi }}_{-1}\), then \(\mu _0(k+{\bar{p}}^\omega )=\mu (k+{\bar{p}}^\omega )-\mu ({\bar{p}}^\omega )+a_\omega k_1 +b_\omega k_2\), see (6.3). Moreover, by using (6.4)–(6.5), we can also rewrite \(\mu (k+{\bar{p}}^\omega )-\mu ({\bar{p}}^\omega )+a_\omega k_1 +b_\omega k_2=\bar{D}_\omega (k)+r_\omega (k)\), which implies that (6.34) at \(h=-1\) is the same as (6.18).

To see that (6.35) holds for \(h=-1\), note that it is equivalent to requiring that \(V^{(-1)}\) is invariant under the transformation \(\varphi ^\pm _{\omega ,x}\rightarrow \varphi ^\pm _{-\omega ,x}\) together with complex conjugation of the kernels. On the other hand, by Remark 6, we know that the potential \(V^{(0)}(\psi ,J)\) is invariant under conjugation of the kernels together with the transformation \(\psi ^\pm _x=(\psi _x^{(0)\pm }+\sum _\omega e^{\pm i {\bar{p}}^\omega }\varphi ^\pm _{\omega ,x})\rightarrow (-1)^x\psi ^\pm _x\), i.e., \(\psi _x^{(0)\pm }\rightarrow (-1)^x\psi _x^{(0)\pm } ,\varphi ^\pm _{\omega ,x}\rightarrow \varphi ^\pm _{-\omega ,x}\). The statement (6.35) for \(h=-1\) easily follows from the relation (6.23) between \(V^{(0)}\) and \(V^{(-1)}\) together with the fact that the propagator \(g^{(0)}\) in (6.17) satisfies

$$\begin{aligned}{}[g^{(0)}(x,y)]^*=(-1)^{x+y} g^{(0)}(x,y), \end{aligned}$$

because \({\bar{p}}^++{\bar{p}}^-=(\pi ,\pi )\).

6.2.3 The inductive step

We assume that (6.32) holds with \(V^{(h)}\) satisfying the properties specified in the inductive statement, and we discuss here how to get the same representation at the next scale \(h-1\). First, we split \(V^{(h)}\) into its local and irrelevant parts: \(V^{(h)}={{\mathcal {L}}}V^{(h)}+{{\mathcal {R}}}V^{(h)}\) where, denoting by \({{\hat{W}}}^{(h),\infty }_{n,m;{\underline{\omega }},{\underline{r}}}\) the infinite volume limit of \({{\hat{W}}}^{(h)}_{n,m;{\underline{\omega }},{\underline{r}}}\),

$$\begin{aligned}&{{\mathcal {L}}}V^{(h)}(\varphi ,J)\nonumber \\&\qquad :=L^{-2}\sum _{\omega }\sum _{k\in {\mathcal {P}}_\omega ({\varvec{\theta }})} {\hat{\varphi }}^+_{k,\omega } [{{\hat{W}}}^{(h),\infty }_{2,0;(\omega ,\omega )}(0)+k\cdot \partial _{k}{{\hat{W}}}^{(h),\infty }_{2,0;(\omega ,\omega )}({0})\big ]{\hat{\varphi }}^-_{k,\omega } \nonumber \\&\quad \qquad +\sum _{x\in \Lambda }\sum _{\omega _1,\ldots ,\omega _4} \varphi ^+_{x,\omega _1}\varphi ^-_{x,\omega _2}\varphi ^+_{x,\omega _3}\varphi ^-_{x,\omega _4} {{\hat{W}}}^{(h),\infty }_{4,0;(\omega _1,\ldots ,\omega _4)}(0,0,0)\nonumber \\&\quad \qquad + \sum _{x\in \Lambda }\sum _{\omega _1,\omega _2,r} J_{x,r}\varphi ^+_{x,\omega _1}\varphi ^-_{x,\omega _2}e^{i({\bar{p}}^{\omega _1} -{\bar{p}}^{\omega _2})x}{{\hat{W}}}^{(h),\infty }_{2,1;(\omega _1,\omega _2),r} (0,{\bar{p}}^{\omega _1}-{\bar{p}}^{\omega _2}).\quad \quad \end{aligned}$$
(6.37)

Remark 10

A few remarks about this definition are in order:

  1. (1)

    The existence of the limit of \({{\hat{W}}}^{(h)}_{n,m;{\underline{\omega }},{\underline{r}}}\) as \(L\rightarrow \infty \) is a corollary of the inductive bounds on the kernels of \(V^{(h)}\), which are uniform in L, as it was the case for \(h=-1\), cf. with Remark 8. More details on the inductive bounds on the kernels of \(V^{(h)}\) are discussed below.

  2. (2)

    The reason why, in the second line of (6.37), we only include terms where the Grassmann fields have the same index \(\omega \), is that the terms with opposite \(\omega \) indices give zero contribution to the generating function, due to the support properties of the Grassmann fields. In fact, in (6.32) we need to compute \(V^{(h)}\) at Grassmann fields \({\hat{\psi }}^{(\le h)\pm }_{k,\omega }\) that, in momentum space, have the same support as \({{\hat{g}}}^{(\le h)}_\omega (k)\), i.e., \(|{\mathcal {M}}^{-1}k|\le c_02^{h}\) (note that the support properties of \({{\hat{g}}}^{(\le h)}_\omega \) are the same as those of \({\bar{\chi }}_h\) (cf. (6.6)), and these were discussed in the third item after (6.3)). If \(h\le -1\) and \(c_0\) is sufficiently small, quadratic terms of the form \({\hat{\psi }}^{(\le h)+}_{k,\omega }{\hat{\psi }}^{(\le h)-}_{k+{\bar{p}}^\omega -{\bar{p}}^{-\omega },-\omega }\) would involve two fields that cannot both satisfy this support property.

  3. (3)

    Due to the Grassmann anti-commmutation rules and the anti-symmetry of the kernels, the quartic term in (6.37) can be rewritten as

    $$\begin{aligned} 4\sum _{x\in \Lambda }\varphi ^+_{x,+}\varphi ^-_{x,+}\varphi ^+_{x,-}\varphi ^-_{x,-} {{\hat{W}}}^{(h),\infty }_{4,0;(+,+,-,-)}(0,0,0). \end{aligned}$$
    (6.38)

Along the induction step, we will need a function \(W^{(h),R}_{2,0;(\omega ,\omega )}(x_1,x_2)\) (the upper index ‘R’ stands for “relativistic”) which should be thought of as the kernel for \(n=2,m=0\) of a relativistic model. More precisely, at step \(h=-1\), one simply lets \(W^{(-1),R}_{2,0;(\omega ,\omega )}(x_1,x_2)\equiv 0\). For \(h<-1\), \(W^{(h),R}_{2,0;(\omega ,\omega )}\) is defined as a suitable modification of \(W^{(h),\infty }_{2,0;(\omega ,\omega )}\) (that, by the induction hypothesis, has already been defined); more precisely, \(W^{(h),R}_{2,0;(\omega ,\omega )}\) is obtained by making the following replacements in \(W^{(h),\infty }_{2,0;(\omega ,\omega )}\) (which should be thought of as a function of the running coupling constants, of the single scale propagators and of the irrelevant part of \(V^{(-1)}\), as explained in Remark 9):

  1. (i)

    the running coupling constants \(\{\nu _{h',\omega },a_{h',\omega },b_{h',\omega }\}_{h'> h}\) are set zero, (note that the running coupling constants \(\lambda _{h'}\) are not set equal to zero);

  2. (ii)

    the single-scale propagators \(g^{(h')}_\omega /Z_{h'-1}\) are replaced by the ‘relativistic’ single-scale propagators \(g^{(h')}_{R,\omega }/Z_{h'-1}\), for all \(h< h'\le -1\), where

    $$\begin{aligned} g^{(h')}_{R,\omega }(x,y)=\int _{{\mathbb {R}}^2} \frac{dk}{(2\pi )^2}e^{-ik(x-y)}\frac{f_{h'}(k)}{{\bar{D}}_\omega (k)}; \end{aligned}$$
    (6.39)
  3. (iii)

    \({\mathcal {R}} V^{(-1)}\) is set to zero.

The function \(W^{(h),R}_{2,0;(\omega ,\omega )}\) will be shown to satisfy both the identity (6.35) and the extra symmetries (in Fourier space)

$$\begin{aligned}&{{\hat{W}}}^{(h),R}_{2,0;(-\omega ,-\omega )}(k)=-[{{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(k)]^*,\nonumber \\&{{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(A^{-1}\sigma _1 Ak)=i\omega [{{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(k)]^*,\nonumber \\&{{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(A^{-1}\sigma _3 Ak)=[{{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(k)]^* \end{aligned}$$
(6.40)

where \(A=\begin{pmatrix} {\bar{\alpha }}^1 &{} {\bar{\beta }}^1 \\ {\bar{\alpha }}^2 &{} {\bar{\beta }}^2 \end{pmatrix}\) while \(\sigma _1,\sigma _3\) are the first and third Pauli matrices. Let us assume that \(W^{(h'),R}_{2,0;(\omega ,\omega )},h'\ge h\) has been already shown to satisfy (6.40) and below we explain how to prove the same at scale \(h-1\).

In order to define the running coupling constants on scale h, we decompose the term containing \( \partial _{k}{{\hat{W}}}^{(h),\infty }_{2,0;(\omega ,\omega )}({0})\) in (6.37), by rewriting

$$\begin{aligned} \partial _{k}{{\hat{W}}}^{(h),\infty }_{2,0;(\omega ,\omega )}({0}) =\partial _{k}{{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}({0}) +\partial _{k}{{\hat{W}}}^{(h),s}_{2,0;(\omega ,\omega )}({0}), \end{aligned}$$
(6.41)

(’s’ stands for ‘subdominant’). From the symmetries (6.40), a straightforward computation (see Appendix 6.6) shows that

$$\begin{aligned} k\cdot \partial _{k}{{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}({0})=-z_h({\bar{\alpha }}_\omega k_1+{\bar{\beta }}_\omega k_2)= -z_h {\bar{D}}_\omega (k), \end{aligned}$$
(6.42)

for some real number \(z_h\). We now combine this term with the Grassmann Gaussian integration \(P_{(\le h)}(d\psi )\), and define:

$$\begin{aligned} P_{(\le h)}(d\psi )e^{-z_hZ_h L^{-2}\sum _\omega \sum _{k\in {\mathcal {P}}'_\omega ({\varvec{\theta }})} {\bar{D}}_\omega (k){\hat{\psi }}^+_{k,\omega }{\hat{\psi }}^-_{k,\omega }}\equiv e^{L^2 t_h}{{\tilde{P}}}_{(\le h)}(d\psi ), \end{aligned}$$
(6.43)

where \({{\tilde{P}}}_{(\le h)}(d\psi )\) is the Grassmann Gaussian integration with propagator

$$\begin{aligned} \frac{{{\tilde{g}}}^{(\le h)}_\omega (x,y)}{Z_{h-1}}=\frac{1}{L^2}\sum _{k\in {\mathcal {P}}_\omega '({\varvec{\theta }})}e^{-ik(x-y)}\frac{{\bar{\chi }}_{h}(k)}{{{\tilde{Z}}}_{h-1}(k) {\bar{D}}_\omega (k)+r_\omega (k)}, \end{aligned}$$
(6.44)

with

$$\begin{aligned} {{\tilde{Z}}}_{h-1}(k):=Z_h(1+z_h{\bar{\chi }}_h(k)),\qquad Z_{h-1}:={{\tilde{Z}}}_{h-1}(0)=Z_h(1+z_h), \end{aligned}$$
(6.45)

and \(e^{L^2 t_h}\) is a constant that normalizes \({{\tilde{P}}}_{(\le h)}(d\psi )\) to 1:

$$\begin{aligned} t_h=\frac{1}{L^2}\sum _\omega \sum _{k\in {\mathcal {P}}'_\omega ({\varvec{\theta }})} \log \Big (1+\frac{z_h{\bar{\chi }}_h(k){\bar{D}}_\omega (k)}{{\bar{D}}_\omega (k)+r_\omega (k)/Z_h}\Big ) . \end{aligned}$$
(6.46)

By using (6.43), we rewrite the Grassmann integral in the right side of (6.32) as

$$\begin{aligned}&\int P_{(\le h)}(d\psi )e^{V^{(h)}(\sqrt{Z_h}(\psi +\Psi ),J)}=e^{L^2 t_h-z_hZ_h L^{-2}\sum _\omega {\bar{D}}_\omega (k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega ){\hat{\Psi }}^+_{\omega } {\hat{\Psi }}^-_{\omega }}\nonumber \\&\quad \times \int {{\tilde{P}}}_{(\le h)} (d\psi )e^{{{\widehat{V}}}^{(h)}(\sqrt{Z_{h-1}} (\psi +\Psi ),J)}, \end{aligned}$$
(6.47)

where

$$\begin{aligned} {{\widehat{V}}}^{(h)}(\varphi ,J)= & {} L^{-2}\sum _{\omega }\sum _{k\in {\mathcal {P}}_\omega ({\varvec{\theta }})} {\hat{\varphi }}^+_{k,\omega } [2^h\nu _{h,\omega }+a_{h,\omega }k_1+b_{h,\omega }k_2\big ]{\hat{\varphi }}^-_{k,\omega } \nonumber \\&+\lambda _h \sum _{x\in \Lambda } \varphi ^+_{x,+}\varphi ^-_{x,+}\varphi ^+_{x,-}\varphi ^-_{x,-} \nonumber \\&+ \sum _{\omega _1,\omega _2,r}\frac{Y_{h,r,(\omega _1,\omega _2)}}{Z_{h-1}}\sum _{x\in \Lambda } J_{x,r}e^{i({\bar{p}}^{\omega _1}-{\bar{p}}^{\omega _2})x} \varphi ^+_{x,\omega _1}\varphi ^-_{x,\omega _2}\nonumber \\&+{\mathcal {R}}V^{(h)}(\sqrt{{Z_h}/{Z_{h-1}}}\,\varphi ,J), \end{aligned}$$
(6.48)

and the running coupling constants at scale h are defined as

$$\begin{aligned}&2^h\nu _{h,\omega }=\frac{Z_h}{Z_{h-1}} {{\hat{W}}}^{(h),\infty }_{2,0;(\omega ,\omega )}(0),\nonumber \\&a_{h,\omega }=\frac{Z_h}{Z_{h-1}} \partial _{k_1}{{\hat{W}}}^{(h),s}_{2,0;(\omega ,\omega )}(0), \qquad b_{h,\omega }=\frac{Z_h}{Z_{h-1}} \partial _{k_2}{{\hat{W}}}^{(h),s}_{2,0;(\omega ,\omega )}(0), \nonumber \\&\lambda _h =4\Big (\frac{Z_h}{Z_{h-1}}\Big )^2 {{\hat{W}}}^{(h),\infty }_{4,0;(+,+,-,-)}(0,0,0),\nonumber \\&Y_{h,r,(\omega _1,\omega _2)}=Z_h {{\hat{W}}}^{(h),\infty }_{2,1;(\omega _1,\omega _2),r} (0,{\bar{p}}^{\omega _1}-{\bar{p}}^{\omega _2}). \end{aligned}$$
(6.49)

Thanks to the symmetry (6.35) of the kernels (that by inductive hypothesis holds at step h) the running coupling constants satisfy the following:

$$\begin{aligned} \nu _{h,\omega }=\nu _{h,-\omega }^*, \quad a_{h,\omega }=-a_{h,-\omega }^*, \quad b_{h,\omega }=-b_{h,-\omega }^*, \quad Y_{h,r,{\underline{\omega }}}=Y_{h,r,-{\underline{\omega }}}^*. \end{aligned}$$
(6.50)

Moreover \(\lambda _h\in {\mathbb {R}}\): for this, one uses both (6.35) and the fact that

$$\begin{aligned} {{\hat{W}}}^{(h)}_{4,0;(+,+,-,-)}(0,0,0)={{\hat{W}}}^{(h)}_{4,0;(-,-,+,+)}(0,0,0). \end{aligned}$$

For later reference, we rewrite the local part of \({{\widehat{V}}}^{(h)}(\varphi ,J)\) as

$$\begin{aligned} {\mathcal {L}}{{{\widehat{V}}}}^{(h)}(\varphi ,J)= & {} \sum _{\omega }\Big [2^h\nu _{h,\omega }F_{\nu ;\omega }(\varphi )+a_{h,\omega }F_{a;\omega }(\varphi )+b_{h,\omega } F_{b;\omega }(\varphi )\Big ]\nonumber \\&+\lambda _h F_\lambda (\varphi )+ \sum _{r,{\underline{\omega }}}\frac{Y_{h,r,{\underline{\omega }}}}{Z_{h-1}}F_{Y;r,{\underline{\omega }}}(\varphi ,J), \end{aligned}$$
(6.51)

(for the definitions of \(F_{\nu ;\omega }(\varphi ),F_{a;\omega }, F_{b;\omega }\), etc., compare (6.51) with the first two lines of (6.48)).

We now decompose the propagator (6.44) as

$$\begin{aligned} {{\tilde{g}}}^{(\le h)}_\omega (x,y)= g^{(h)}_\omega (x,y)+g^{(\le h-1)}_\omega (x,y), \end{aligned}$$

with \(g^{(\le h-1)}_\omega \) as in (6.34) and \(g^{(h)}_\omega \) as in (6.36). To see that this decomposition holds, note that \({{\tilde{Z}}}_{h-1}(k)\equiv Z_{h-1}\) on the support of \({\bar{\chi }}_{h-1}(\cdot )\).

Then, rewrite (6.47) as

$$\begin{aligned}&\int P_{(\le h)}(d\psi )e^{V^{(h)}(\sqrt{Z_h}(\psi +\Psi ),J)}=e^{L^2 t_h-z_hZ_h L^{-2}\sum _\omega {\bar{D}}_\omega (k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega ) {\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_{\omega }}\nonumber \\&\qquad \times \int P_{(\le h-1)}(d\psi ) \int P_{(h)}(d\psi ')e^{{{\widehat{V}}}^{(h)}(\sqrt{Z_{h-1}} (\psi +\psi '+\Psi ),J)}, \end{aligned}$$
(6.52)

which implies the validity of the representation (6.32) at scale \(h-1\), with \(E^{(h-1)}\), \(S^{(h-1)}(\cdot )\) and \(V^{(h-1)}(\cdot )\) defined by

$$\begin{aligned}&e^{L^2 E^{(h-1)} + S^{(h-1)}(J)+ V^{(h-1)}(\sqrt{Z_{h-1}}(\psi +\Psi ),J)}\nonumber \\&\qquad =e^{L^2 (E^{(h)}+t_h)+S^{(h)}(J)}\int P_{(h)}(d\psi ') e^{{{\widehat{V}}}^{(h)}(\sqrt{Z_{h-1}}(\psi +\psi '+\Psi ),J)}, \end{aligned}$$
(6.53)

that is,

$$\begin{aligned}&L^2 (E^{(h-1)}-E^{(h)}-t_h) +(S^{(h-1)}(J)-S^{(h)}(J))+ V^{(h-1)}(\varphi ,J) \nonumber \\&\quad =\sum _{n\ge 1}\frac{1}{n!}{\mathcal {E}}^T_h(\underbrace{{{\widehat{V}}}^{(h)} (\varphi +\sqrt{Z_{h-1}}\psi ',J);\cdots ; {{\widehat{V}}}^{(h)}(\varphi +\sqrt{Z_{h-1}}\psi ',J)}_{n\ \mathrm{times}}, \end{aligned}$$
(6.54)

with \({\mathcal {E}}^T_h\) the truncated expectation w.r.t. the Grassmann Gaussian integration \(P_{(h)}(d\psi )\), and \(E^{(h-1)}\), \(S^{(h-1)}(\cdot )\) fixed as usual by the conditions \(S^{(h-1)}(0)=0\) and \(V^{(h-1)}(0,J)=0\).

To conclude the proof of the induction step, it remains to prove that the kernels of \(V^{(h-1)}\) satisfy (6.35) and that (6.40) holds, at scale \(h-1\). The proof of the former statement is very similar (but not identical) to the argument used in Sect. 6.2.2 to prove (6.35) at scale \(h=-1\) starting from the symmetries of \(V^{(0)}\). Namely, thanks to (6.35) at scale h, the potential \(V^{(h)}\) is invariant under \(\varphi ^\pm _{x,\omega }\rightarrow \varphi ^\pm _{x,-\omega }\) together with complex conjugation of the kernels. Then, the claim follows from the representation (6.54), together with the fact that the propagator \(g^{(h)}\) (defined in (6.36)) satisfies the symmetry

$$\begin{aligned}{}[g^{(h)}_\omega (x,y)]^*=g^{(h)}_{-\omega }(x,y). \end{aligned}$$
(6.55)

As for (6.40) at scale \(h-1\), the proof uses the symmetries of the relativistic propagator (6.39), together with the fact that \(\lambda _{h'}\) is real. See Appendix 6.6.

Remark 11

Note that, if the function \(z_{h'}\) in (6.45) is sufficiently small for all the scales \(h\le h'\le -1\), say \(|z_{h'}|\le \varepsilon \) uniformly in \(L,h'\), then \(e^{-c \varepsilon |h|}\le Z_h\le e^{c \varepsilon |h|}\). As a consequence, \(g^{(h)}_\omega \) satisfies a bound analogous to (6.21), namely

$$\begin{aligned} |g^{(h)}_\omega (x,y)|\le C_0 2^h e^{-\kappa \sqrt{2^h|x-y|}}. \end{aligned}$$
(6.56)

In fact, note that on the support of \(f_h(\cdot )\) (which is concentrated on \(k:|k|\sim 2^h\)), \({{\tilde{Z}}}_{h-1}(k)/Z_{h-1}=1+O(\varepsilon )\) and recall that \(r_\omega (\cdot )\) is quadratic for small values of its argument, so that \(r_\omega (k)/Z_{h-1}\) is negligible w.r.t. \({\bar{D}}_\omega (k)\). The propagator \(g^{(h)}_{R,\omega }\) satisfies the same estimate as (6.56), while the difference \(g^{(h)}_{\omega }-g^{(h)}_{R,\omega }\) satisfies an estimate that is better by a factor \(2^h\).

6.2.4 The beta function

The iterative integration scheme described above allows us to express the kernels of \(V^{(h)}\) and, in particular, the running coupling constants (RCC) at scale h, as functions of the sequence of RCC on higher scales, \(\{\lambda _{h'},\nu _{h',\omega },a_{h',\omega },b_{h',\omega },Y_{h',r,{\underline{\omega }}}\}_{h<h'\le -1}\), of the single-scale propagators \(\{g^{(h')}_\omega /Z_{h'-1}\}_{h<h'\le -1}\), and of \({\mathcal {R}} V^{(-1)}\). That is, we can rewrite Equation (6.49) in the form

$$\begin{aligned}&\nu _{h-1,\omega }=2\nu _{h,\omega }+B^\nu _{h,\omega }\quad a_{h-1,\omega }=a_{h,\omega }+B^a_{h,\omega },\quad b_{h-1,\omega }=b_{h,\omega }+B^b_{h,\omega },\nonumber \\&\lambda _{h-1}=\lambda _h+B^\lambda _h,\qquad Y_{h-1,r,{\underline{\omega }}}=Y_{h,r,{\underline{\omega }}}+B^Y_{h,r,{\underline{\omega }}}, \end{aligned}$$
(6.57)

where \(B^\#_{h,\cdot }\), \(h\le -1\), is the so-called Beta function. One has to think of \(B^\#_{h,\cdot }\) as a function of the RCC on scales \(h'\) with \(h\le h'\le 0\). Note that the first four equations makes sense also with \(h=0\), in which case they express the relation between \((\nu _{-1,\omega },a_{-1,\omega },b_{-1,\omega },\lambda _{-1})\) and \((\nu _{0,\omega },a_{0,\omega },b_{0,\omega },\lambda _{0})\), see (6.29). Note also that by construction the beta function \(B^\#_{h,\cdot }\) depends on \(Z_{h'}\) only via the combinations \(Z_{h'}/Z_{h'-1}=(1+z_{h'})^{-1}\), with \(h< h'< 0\). For later reference, we rewrite the definition of \(z_h\), (6.42), in a form analogous to (6.57),

$$\begin{aligned} z_{h-1}=B^z_h, \quad h\le 0, \end{aligned}$$
(6.58)

where the right side is thought of as a function of \((\lambda _{h'},z_{h'})_{h\le h'\le 0}\), with the convention that \(z_{0}=z_{-1}=0\) (the latter is because \(W^{(-1),R}_{2,0;(\omega ,\omega )}\equiv 0\)).

Remark 12

The components of the beta function for \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega },\lambda _h\) are independent of \(Y_{h',r,{\underline{\omega }}},h'>h\). Therefore, we can first solve the flow equation for \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega },\lambda _h\) and then inject the solution into the flow equation for \(Y_{h,r,{\underline{\omega }}}\).

Before we proceed in describing the dimensional bounds satisfied by the kernels of the effective potential, let us comment on their structure. We have proven inductively that \(V^{(h)}\) has, in momentum space, the same structure as (6.26). If one writes \(V^{(h)}\) in real space, due to iterative action of the \({\mathcal {R}}\) operator in the inductive procedure explained above, the structure that emerges naturally is that of a polynomial with spatial derivatives acting on some of the Grassmann fields \(\varphi ^\pm _{x,\omega }\). For an explanation of why this is the case see [25, Section 6.1.4] and Appendix B below, where finite-size effects associated with the action of \({\mathcal {R}}\) are also discussed. Correspondingly, \(V^{(h)}\) can be represented as

$$\begin{aligned}&V^{(h)}(\varphi ,J)= \sum _{\begin{array}{c} n,m\ge 0:\\ n\ {\mathrm{even}}, \ n\ge 2 \end{array}}\sum _{{\underline{x}},\, {\underline{y}},\, {\underline{\omega }},\,{\underline{r}},\, {\underline{i}},{\underline{q}}} W^{(h)}_{n,m,{\underline{i}},{\underline{q}};{\underline{\omega }},{\underline{r}}}({\underline{x}}, {\underline{y}})\nonumber \\&\qquad \times \, {\hat{\partial }}^{q_{1}}_{i_1}\varphi ^{(\le h)+}_{x_1,\omega _1}\cdots {\hat{\partial }}^{q_{n}}_{i_n}\varphi ^{(\le h)-}_{x_n,\omega _n} J_{y_1,r_1}\cdots J_{y_m,r_m}. \end{aligned}$$
(6.59)

The main difference between this formula and (6.24), besides the scale label h replacing \(-1\), is the presence of the indices \({\underline{i}}=(i_1,\ldots ,i_n)\in \{1,2\}^n\) and \({\underline{q}}=(q_1,\ldots ,q_n)\in \{0,1,2\}^n\) and the operators \({\hat{\partial }}^{q}_i\) acting on the Grassmann fields: this is a differential operator, dimensionally equivalent to a derivative of order q in direction i. Let us stress that the representation in (6.59) is not unique: the claim is that there exists such a representation, with the kernels satisfying natural dimensional estimates, discussed below.

In order for the iterative construction to allow us to compute the thermodynamic and correlation functions, we need to prove that: (i) the RCC \(\nu _{h,\omega },a_{h,\omega }, b_{h,\omega },\lambda _{h},z_h\) are small, uniformly in the scale (say, smaller than a sufficiently small constant \(\varepsilon \)), provided the functions \({\bar{p}}^\omega ,a_\omega ,b_\omega \) (see (6.3)) have been properly fixed; (ii) the kernels of the effective potential are all well defined (i.e. the sums (6.54) are convergent uniformly in L), quasi-local (i.e., fast decaying, with a stretched-exponential behavior) and satisfy natural scaling properties, i.e.,

$$\begin{aligned} \Vert W^{(h)}_{n,m,{\underline{i}}, {\underline{q}}}\Vert _{\kappa ,h}\le C^{n+m} \varepsilon ^{\max \{1,cn\}} 2^{h(2-n/2-m-|{\underline{q}}|)} \Big (\max _{h'\ge h}\frac{|Y_{h',\cdot }|}{|Z_{h'}|}\Big )^m, \end{aligned}$$
(6.60)

with \(|{\underline{q}}|=\sum _{i=1}^nq_{i}\), \(|Y_{h',\cdot }|=\max _{r,{\underline{\omega }}}|Y_{h',r,{\underline{\omega }}}|\), and

$$\begin{aligned} \Vert W^{(h)}_{n,m,{\underline{i}}, {\underline{q}}}\Vert _{\kappa ,h}:=L^{-2}\sup _{{\underline{\omega }}, {\underline{r}}}\sum _{{\underline{x}}, {\underline{y}}}|W^{(h)}_{n,m,{\underline{i}}, {\underline{q}};{\underline{\omega }},{\underline{r}}}({\underline{x}}, {\underline{y}})|e^{\kappa \sqrt{2^{h}d({\underline{x}},{\underline{y}})}}, \end{aligned}$$
(6.61)

for suitable constants \(C,c,\kappa \), independent of Lh.

The boundedness of the flow of the RCC and the validity of the dimensional bounds for the kernels of the effective potential will, in fact, be the final outcome of our analysis. The logic of proof goes as follows: one first proves the validity of the dimensional bounds on the kernels, under the assumption that the RCC remain small. These bounds will, in particular, imply that the components of the beta function are well defined and satisfy bounds that are uniform in L and h. This part of the proof is pretty standard: it follows from a representation of the effective potential in terms of Gallavotti-Nicolò (GN) trees, see Sect. 6.3 below, and an iterative application of the Battle-Brydges-Federbush-Kennedy (BBFK) determinant formula, see, e.g., [25, Lemma 3].

Next, we prove that the RCC remain bounded, by studying the flow generated by the beta function. The key point is that, as long as the RCC on scales larger than h are small, then the beta function on scale h is well defined, and can be used to control the evolution of the RCC for another step. This opens the way to an inductive proof of the smallness of the RCC. Of course, the fact that RCC remain small at all scale requires a specific (model-dependent) structure of the beta function. In our case, we are lucky enough that the beta function has structure which maintains the RCC small at all scale, provided the initial data are small, and that \({\bar{p}}^\omega ,a_\omega ,b_\omega \) are properly fixed, see Sect. 6.4 below. It is not just a matter of luck, of course: a key point in the analysis is played by the comparison of the \(\lambda \)-component of the beta function of our dimer model, with the corresponding quantity for the reference continuum model (the two functions are the same at dominant order). The exact solvability of the reference model implies the validity of a remarkable cancellation in the \(\lambda \)-component of the beta function for the reference model and, therefore, a posteriori, for our dimer model, as well.

6.3 Tree expansion for the effective potential

As anticipated above, the detailed structure of the kernels of \(V^{(h)}\), arising from the iterative construction described in the previous section, can be conveniently represented in terms of GN trees. The definition of GN trees, of their values, and the procedure leading to their introduction have been discussed at length in several previous papers and will not be repeated here, see e.g. [21]; in particular, we refer to [25, Section 5.2.1 and 6.2] for a description of the GN tree expansion in a context very similar to the present one, i.e., in the case of isotropic, ‘tilt-less’, interacting dimer models with weights \({\underline{t}}\equiv 1\) and plaquette interaction. The present case differs from the one treated in [25] for the fact that here the model is anisotropic (and, correspondingly, the height has an average slope that is different from zero). Technically, this means that in the present case the expansion involves more running coupling constants than those considered in [25]: the RCC \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega }\) are identically zero in the tilt-less case. In particular, the trees involved in our construction are characterized by the following features, slightly different from those listed in [25, Section 6.2]:

  1. (1)

    A GN tree \(\tau \) contributing to \(V^{(h)}\), \({{\tilde{S}}}^{(h)}(J):=S^{(h)}(J)-S^{(h+1)}(J)\), or to \({{\tilde{E}}}^{(h)}=E^{(h)}-E^{(h+1)}-t_{h+1}\) has root on scale h and can have endpoints (either ‘normal’ or ‘special’, which are those represented as black dots or white squares, respectively, in [25], see, e.g., [25, Fig.13]) on all possible scales between \(h+2\) and 0. The endpoints v on scales \(h_v< 0\) are preceded by a node \(v'\) of \(\tau \), on scale \(h_{v'}=h_v-1\), that is necessarily a branching point. The family of GN trees with root on scale h, \(N_n\) normal endpoints and \(N_s\) special endpoints is denoted by \({\mathcal {T}}^{(h)}_{N_n,N_s}\).

  2. (2)

    A normal endpoint v on scale \(h_v\le 0\) can be of five different types, \(\lambda ,\nu ,a,b\), or \({\mathcal {R}} V^{(-1)}\). If v is of type \(\lambda ,\nu ,a\) or b, then it is associated with \(\lambda _{h_{v'}} F_\lambda (\sqrt{Z_{h_{v'}-1}}\varphi ^{(\le h_{v'})})\), or \(\sum _\omega \nu _{h_{v'},\omega }F_{\nu ;\omega }(\sqrt{Z_{h_{v'}-1}}\varphi ^{(\le h_{v'})})\), or \(\sum _\omega a_{h_{v'},\omega }F_{a;\omega }(\sqrt{Z_{h_{v'}-1}}\varphi ^{(\le h_{v'})})\), or \(\sum _\omega b_{h_{v'},\omega }F_{b;\omega }(\sqrt{Z_{h_{v'}-1}}\varphi ^{(\le h_{v'})})\), depending on its type (recall that the monomials \(F_\lambda \), \(F_{\nu ;\omega }\), etc., were defined in (6.51)); in this case, the node \(v'\) immediately preceding v on \(\tau \), of scale \(h_{v'}=h_v-1\), is necessarily a branching point. If v is of type \({\mathcal {R}} V^{(-1)}\), then \(h_v=0\), and v is associated with (one of the monomials contributing to) \(\mathcal RV^{(-1)}(\varphi ^{(\le -1)},0)\); in this case, the node immediately preceding v on \(\tau \), of scale \(h_v-1\), is not necessarily a branching point.

  3. (3)

    A special endpoint v on scale \(h_v\le 0\) can be either local, or non-local. If v is local, then it is associated with

    $$\begin{aligned} \frac{Y_{h_{v'},r,{\underline{\omega }}}}{Z_{h_{v'}-1}}F_{Y;r,{\underline{\omega }}}(\sqrt{Z_{h_{v'}-1}}\varphi ^{(\le h_{v'})},J), \end{aligned}$$
    (6.62)

    for some \(r\in \{1,2,3,4\}\), \({\underline{\omega }}=(\omega _1,\omega _2)\in \{\pm \}^2\); if \(\omega _1=\omega _2\), we shall say that v is a ‘density endpoint’, while, if \(\omega _1\ne \omega _2\), that v is a ‘mass endpoint’. Note that the factors \(Z_{h_{v'}-1}\) in (6.62) simplify: the summand equals \(Y_{h_{v'},r,{\underline{\omega }}}F_{Y;r,{\underline{\omega }}}(\varphi ^{(\le h_{v'})},J)\); in (6.62), these factors are kept just for uniformity of notation with the cases in the previous item. In the case that v is local, the node \(v'\) immediately preceding v on \(\tau \), of scale \(h_{v'}=h_v-1\), is necessarily a branching point. If v is non-local, then \(h_v=0\), and v is associated with (one of the monomials contributing to) \(V^{(-1)}(\varphi ^{(\le -1)},J)-V^{(-1)}(\varphi ^{(\le -1)},0)\); in this case, the node immediately preceding v on \(\tau \), of scale \(h_v-1\), is not necessarily a branching point.

In addition to the items above, let us recall that each vertex of the tree that is not an endpoint and that is not the special vertex \(v_0\) (the leftmost vertex of the tree, immediately following the root on \(\tau \)) is associated with the action of an \({{\mathcal {R}}}\) operator.

In terms of the tree expansion, we can express the effective potential and the single-scale contributions to the free energy and generating function as

$$\begin{aligned} L^2 {{\tilde{E}}}^{(h)}+{{\tilde{S}}}^{(h)}(J)+V^{(h)}(\sqrt{Z_h}\varphi ,J)=\sum _{\begin{array}{c} N_n,N_s\ge 0:\\ N_n+N_s\ge 1 \end{array}}\sum _{\tau \in {\mathcal {T}}^{(h)}_{N_n,N_s}} V^{(h)}(\tau ,\sqrt{Z_h}\varphi ,J), \end{aligned}$$
(6.63)

where

$$\begin{aligned}&V^{(h)}(\tau ,\sqrt{Z_h}\varphi ,{J})=\nonumber \\&\quad =\sum _{ \mathbf{P}\in {{\mathcal {P}}}_\tau }\sqrt{Z_h}^{|P_{v_0}^\psi |}\sum _{T\in \mathbf{T}}\sum _\mathbf{i, \mathbf{q}}\sum _{ \mathbf{x}_{v_0}}W_{\tau ,\mathbf{P},{T},\mathbf{i}, \mathbf{q}} (\mathbf{x}_{v_0}) D_\mathbf{i}^\mathbf{q}\varphi ({P_{v_0}^\psi })\, J({P_{v_0}^J}) \;. \end{aligned}$$
(6.64)

Equation (6.64) is the analogue of [25, (6.64)], and we refer the reader to that paper for the notation and a sketch of the proof (in this formula, the indices \(\mathbf{i},\mathbf{q}\) replace the multi-indices \(\beta \in B_T\) [25, (6.64)]). We recall that \(P_{v_0}^\psi \) and \(P_{v_0}^J\) are two sets of indices that label the Grassmann external fields and the external fields of type J, respectively; moreover, \(J({P_{v_0}^J})=\prod _{f\in P_{v_0}^J}J_{y(f),r(f)}\) and

$$\begin{aligned} D_\mathbf{i}^\mathbf{q}\varphi ({P_{v_0}^\psi })=\prod _{f\in P_{v_0}^\psi } {\hat{\partial }}_{i(f)}^{q(f)}\varphi _{x(f),\omega (f)}^{\varepsilon (f)}\;. \end{aligned}$$
(6.65)

Clearly, the kernels in (6.59) are obtained by summing \(W_{\tau ,\mathbf{P},{T},\mathbf{i},\mathbf{q}} (\mathbf{x}_{v_0})\) over \(\tau \in {\mathcal {T}}^{(h)}_{N_n,N_s}\) and over \(N_n,N_s\), under the constraint that the number of external fields of type \(\psi \) and J is equal to n and m, respectively, that the elements of \(\mathbf{i}\) are the same as \({\underline{i}}\), etc. Similarly, the kernels of the single scale contribution to the generating function, \({{\tilde{S}}}^{(h)}(J)\), which we denote by \(W^{(h)}_{0,m;{\underline{r}}}({\underline{y}})\), are obtained by summing the tree values \(W_{\tau ,\mathbf{P},{T}}\) over \(\tau \in {\mathcal {T}}^{(h)}_{N_n,N_s}\) and over \(N_n,N_s\), under the constraint that \(P_{v_0}^\psi =\emptyset \) and that \(\cup _{f\in P_{v_0}^J}\{(y(f), r(f))\}\) matches the tuple \(({\underline{y}},{\underline{r}})\); finally, the single scale contribution to the free energy, \(L^2 {{\tilde{E}}}^{(h)}\) is obtained by an analogous sum over GN trees, under the constraint that \(P_{v_0}^\psi =P_{v_0}^J=\emptyset \).

The bound (6.60), as well as the analogous one for \(W_{0,m:{\underline{r}}}\), is a corollary of the following fundamental bound on the weighted \(L_1\) norm of \(W_{\tau ,\mathbf{P},{T},\mathbf{i},\mathbf{q}}\), which is the analogue of [25, Proposition 8] and of [11, (3.110)]; for the proof, we refer the reader to [11, 25]. See also Appendix  B below for some technical details.

Proposition 3

There exists L-independent constants \(\varepsilon ,C,c,\kappa >0\) such that, if

$$\begin{aligned} \max _{h'>h}\{|\lambda _{h'}|,|\nu _{h',\omega }|,|a_{h',\omega }|,|b_{h',\omega }|,|z_{h'}|\}\le \varepsilon , \end{aligned}$$
(6.66)

and \(\tau \in {\mathcal {T}}^{(h)}_{N_n,N_s}\), then

$$\begin{aligned}&\Vert W_{\tau ,\mathbf{P},T,\mathbf{i},\mathbf{q}}\Vert _{\kappa ,h}\le C^{N_s} \,(C\varepsilon )^{\max \{N_n,c|I_{v_0}^\psi |\}}\, 2^{h(2-\frac{1}{2}|P_{v_0}^\psi |- |P_{v_0}^J|-|\mathbf{q}|)}\nonumber \\&\quad \times \Biggl [\prod _{v\ {\mathrm{s.e.p.}}}\sup _{r,{\underline{\omega }}} \Big |\frac{Y_{h_v-1,r,{\underline{\omega }}}}{Z_{h_v-1}}\Big |\Biggr ]\Big [\prod _{\begin{array}{c} v\,{\mathrm{not}}\\ {{\mathrm{e.p.}}} \end{array}}\quad \frac{C^{\sum _{i=1}^{s_v}|P_{v_i}|-|P_v|}}{s_v!} \;2^{ \frac{\varepsilon }{2}|P^\psi _v|}2^{2-\frac{1}{2}|P_v^\psi |-|P_v^J|-z(P_v)} \Big ]\;,\nonumber \\ \end{aligned}$$
(6.67)

where: \(|I_{v_0}^\psi |=\sum _{v\ {\mathrm{e.p.}}}|P_{v}^\psi |\) is the total number of Grassmann fields associated with the endpoints of the tree; the first product in the second line runs over the special endpoints, while the second over all the vertices of the tree that are not endpoints. Moreover \(|\mathbf{q}|=\sum _{f\in P_{v_0}^\psi }q(f)\) and

(6.68)

The dimensional gain \(2^{-z(P_v)}\) associated with the marginal and relevant vertices, i.e., those with \((|P_v^\psi |,|P_v^J|)=(4,0),(2,1),(2,0)\), comes from the action of \({\mathcal {R}}\), as explained in [25, Section 6.1.4] and in Appendix B below.

Since the exponents \(2+ \frac{\varepsilon }{2}|P^\psi _v|-\frac{1}{2}|P_v^\psi |-|P_v^J|-z(P_v)\) in (6.67) are all strictly negative, one can sum (6.67) over \(\tau \in {\mathcal {T}}^{(h)}_{N_n,N_s}\), over \(T\in \mathbf{T}\), and over \(\mathbf{P}\in {{\mathcal {P}}}_\tau \), under the constraint that \(|P_{v_0}^\psi |=n\) and \(|P_{v_0}^J|=m\), we get the bound (6.60); see also the discussion after [25, Proposition 8]. Similarly, if we sum (6.67) over \(\tau \in {\mathcal {T}}^{(h)}_{N_n,N_s}\), \(T\in \mathbf{T}\), \(\mathbf{P}\in {{\mathcal {P}}}_\tau \), with \(|P_{v_0}^\psi |=n\), \(|P_{v_0}^J|=m\), under the additional constraint that\(\tau \)has at least one node on scale\(k>h\), then we get a bound that is the same as (6.60) times an additional gain factor \(2^{\theta ' (h-k)}\), where \(\theta '\) is a positive constant, smaller than 1 (estimates are not uniform as \(\theta '\rightarrow 1^-\); from here on, we will choose \(\theta '=3/4\)). This is the so-called short memory property, see Remark 16 after [25, Proposition 8].

6.4 The flow of the running coupling constants

As explained above, as long as the RCC \(\nu _{h',\omega },a_{h',\omega }, b_{h',\omega },\lambda _{h'},z_{h'}\) stay small, for all \(h'>h\), in the sense of (6.66), the beta function controlling the flow of the same constants on all scales larger or equal to h, see (6.57)–(6.58), can be represented in terms of a convergent GN expansion, induced by the one of the kernels of the effective potential discussed above. The goal is then to fix the initial data \(\nu _{0,\omega }, a_{0,\omega }, b_{0,\omega }\), in such a way that the resulting flow of \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega },\lambda _h,z_h\) driven by the beta function stays uniformly small in the scale index. For this purpose, not only we have to make a careful choice of the ‘counter-terms’ \(\nu _{0,\omega }, a_{0,\omega }, b_{0,\omega }\), but we also need to exploit a number of remarkable cancellations, some of which follow from the exact solution of the reference model of Sect. 4. Let us emphasize that we have the right to fix the counter-terms, which up to now were chosen arbitrarily in (6.3) (recall (6.29)), but we cannot change \(\lambda _0=\lambda \), which enters the definition of the model.

We look for a solution of the flow equation for the RCC such that, as \(h\rightarrow -\infty \):

  1. (1)

    \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega }\) tend exponentially to zero; more precisely, recalling that \(|\nu _{h,+}|=|\nu _{h,-}|\), and similarly for \(|a_{h,\omega }|\), \(|b_{h,\omega }|\), we require that

    $$\begin{aligned} \Vert ({\underline{\nu }},{\underline{a}},{\underline{b}})\Vert _\theta :=\sup _{h\le 0}\{2^{-\theta h}|\nu _{h,+}|,\, 2^{-\theta h}|a_{h,+}|,\, 2^{-\theta h}|b_{h,+}|\}\le \varepsilon , \end{aligned}$$
    (6.69)

    for \(\varepsilon \) small enough, \(\theta =1/2\) and \(\underline{\nu }:=(\nu _{h,\omega })^{\omega \in \{\pm \}}_{h\le 0}, \underline{a}:=(a_{h,\omega })^{\omega \in \{\pm \}}_{h\le 0}\) and \(\underline{b}:=(b_{h,\omega })^{\omega \in \{\pm \}}_{h\le 0}\).

  2. (2)

    \(\lambda _h\) tends exponentially to a finite limiting value \(\lambda _{-\infty }\); more precisely, given \(\underline{\lambda }:=(\lambda _h)_{h\le -1}\) and a positive constant \(\varepsilon '\) smaller than the constant \(\varepsilon \) in the previous item, we require that \(|\lambda _0|\le \varepsilon '\) and

    $$\begin{aligned} \Vert {\underline{\lambda }}\Vert _\theta :=\sup _{h\le 0}\{2^{-\theta h}|\lambda _{h-1}-\lambda _h|\}\le \varepsilon ', \end{aligned}$$
    (6.70)

    where \(\theta \) is the same as in the previous item; note that, from the condition on \(\lambda _0\) and (6.70), it follows that

    $$\begin{aligned} |\lambda _h|\le \frac{\varepsilon '}{1-2^{-\theta }}+\varepsilon ', \end{aligned}$$
    (6.71)

    uniformly in h.

In order to construct such solution, we proceed in various steps.

  • First, in Sect. 6.4.1, given any sequence \({\underline{\lambda }}\) satisfying the conditions of item (2), we show how to construct the solution \(z_h\) to the beta function equation (6.58) with initial datum \(z_0=z_{-1}=0\).

  • Next, in Sect. 6.4.2, we reformulate the problem of constructing a bounded solution to the beta function equations for \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega },\lambda _h\) into a fixed point problem of the form \({\underline{u}}=T({\underline{u}})\) in a Banach space of sequences.

  • In order to construct and prove uniqueness of such a fixed point we appeal to the Banach fixed point theorem, a.k.a. the contraction theorem: in Sect. 6.4.3 we prove the existence of an invariant set \(X_\varepsilon \), i.e., a set such that \(T(X_\varepsilon )\subseteq X_\varepsilon \), and in Sect. 6.4.4 we prove that T is a contraction on \(X_\varepsilon \).

  • In the remaining subsections, we collect various consequences of the proof of existence of the fixed point (or, equivalently, of the desired bounded sequence of RCCs): in Sect. 6.4.5 we invert the counterterms and explain how to compute the dressed Fermi points \({\bar{p}}^{\omega }\) and Fermi velocities \({\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \) (the functions \({\bar{p}}^{\omega }\), \({\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \) constructed here are the same as items (1)–(2) of Proposition 2); in Sect. 6.4.6 we discuss the flow of \(Z_h\) and define the associated critical exponent \(\eta \); in Sects. 6.4.7 and  6.4.8 we fix the bare couplings \(\lambda _\infty \) and Z of the reference model, in such a way that the flows of the reference and dimer model have the same asymptotic behavior as \(h\rightarrow -\infty \), including the same critical exponents (the choices of \(\lambda _\infty , Z\) made there are those stated in item (3) of Proposition 2); finally, in Sect. 6.4.9, we discuss the flow of \(Y_{h,r,{\underline{\omega }}}\) and define the associated critical exponents \(\eta ,\eta _1\).

6.4.1 Fixing \((z_h)_{h\le -1}\)

Given a sequence \({\underline{\lambda }}:=(\lambda _h)_{h\le -1}\) satisfying (6.70), we construct the solution of the beta function equation

$$\begin{aligned} z_{h-1}=B^z_h({\underline{\lambda }},{\underline{z}}) \end{aligned}$$
(6.72)

for \({\underline{z}}:=(z_h)_{h\le -1}\) iteratively in h, starting from \(h=0\). We denote this solution by \({\underline{z}}^*({\underline{\lambda }})\). We recall that the definition of the function \(B^z_h\), first introduced in (6.58), is induced by the iterative definition of \(z_h\), (6.42). Note that, of course, \(B^z_h({\underline{\lambda }},{\underline{z}})\) only depends on the components of \({\underline{\lambda }},{\underline{z}}\) of scale index larger or equal to h; note also that \(B^z_h\) does not depend on the irrelevant part of the effective interaction and, in particular, it does not depend on \(\lambda _0\), because its definition only involves \(W^{(h),R}_{2,0;(\omega ,\omega )}\) (see (6.42) and the definition of \(W^{(h),R}_{2,0;(\omega ,\omega )}\) given after Remark 10), specifically item (iii) after (6.39)).

By using the tree expansion of the beta function, we now show that the solution \({\underline{z}}^*({\underline{\lambda }})\) of (6.72) is a Cauchy sequence, differentiable in \({\underline{\lambda }}\); more precisely, we prove that, for \(\lambda _0\) fixed, such that \(|\lambda _0|\le \varepsilon '\), and \({\underline{\lambda }}\) satisfying (6.70),

$$\begin{aligned} |z_{h-1}^*({\underline{\lambda }})-z_h^*({\underline{\lambda }})|\le C_0(\varepsilon ')^2 2^{\theta h}, \qquad \Big |\frac{\partial z^*_h ({\underline{\lambda }})}{\partial \lambda _k}\Big |\le C_0\varepsilon ' 2^{\theta (h-k)}, \end{aligned}$$
(6.73)

for all \(h\le k\le -1\). Once this is done, we plug \({\underline{z}}^*({\underline{\lambda }})\) in the flow equations for \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega },\lambda _h\), i.e., the first four equations of (6.57), so that a posteriori their beta functions are re-expressed purely in terms of \(\lambda _0\) and \({\underline{u}}\), where

$$\begin{aligned} {\underline{u}}=({\underline{\nu }}, {\underline{a}}, {\underline{b}},{\underline{\lambda }}), \end{aligned}$$
(6.74)

with \({\underline{\nu }}:=(\nu _{h,\omega })_{h\le 0}^{\omega \in \{\pm \}}\), \({\underline{a}}:=(a_{h,\omega })_{h\le 0}^{\omega \in \{\pm \}}\) and \({\underline{b}}:=(b_{h,\omega })_{h\le 0}^{\omega \in \{\pm \}}\).

Let us prove the first inequality in (6.73), inductively in h. Note that at the first step, \(h=0\), the inequality is trivially true, simply because \(z^*_0({\underline{\lambda }})=z^*_{-1}({\underline{\lambda }})=0\). We assume that \(|z_{h'-1}^*({\underline{\lambda }})-z_{h'}^*({\underline{\lambda }})|\le C_0(\varepsilon ')^2 2^{\theta h'}\), for all scales \(h<h'\le 0\), and we want to prove that the same bound holds for \(h'=h\). Note that, for \(\varepsilon '\) sufficiently small, the first inequality in (6.73) also implies that \(|z^*_{h'}({\underline{\lambda }})|\le \varepsilon '\le \varepsilon \), \(\forall h\le h'\le -1\), uniformly in \({\underline{\lambda }}\). Recall that the definition of \(B^z_{h}\) is induced by (6.42). The kernel \(W^{(h),R}_{2,0;(\omega ,\omega )}\) can be written as a sum over GN trees, analogous to (6.64), and the contribution associated with each tree can be bounded as in Proposition 3. Note that the trees contributing to \(B^z_h\) have only endpoints of type \(\lambda \): this is because \(z_h\) is defined in terms of \(\partial _k {{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(0)\), see (6.42), and \(W^{(h),R}_{2,0;(\omega ,\omega )}\) is obtained, by definition, by: setting the RCCs \(\nu _{h',\omega },a_{h',\omega },b_{h',\omega }\), as well as the irrelevant coupling \({\mathcal {R}} V^{(-1)}\), to zero; replacing the single-scale propagators \(g^{(h')}_{\omega }\) by their ‘relativistic’ counterpart, \(g^{(h')}_{R,\omega }\); see the discussion after (6.38). Moreover, the contribution from the tree with exactly one \(\lambda \) endpoint (which corresponds, in the language of Feynman diagrams, to the ‘tadpole’) is exactly zero (because the kernel of the tadpole is proportional to a delta function, so that the associated contribution to \(\partial _k {{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(0)\) is zero). These considerations lead to the following representation:

$$\begin{aligned} B^z_{h}({\underline{\lambda }},{\underline{z}}^*)=\sum _{N\ge 2}\sum _{\tau \in {\mathcal {T}}^{(h)}_{N,0}} \sum _{\mathbf{P}\in {{\mathcal {P}}}_\tau }\sum _{T\in \mathbf{T}}B^z({\underline{\lambda }},{\underline{z}}^*;\tau ,\mathbf{P},T), \end{aligned}$$
(6.75)

where \({\underline{z}}^*={\underline{z}}^*({\underline{\lambda }})\) (recall that \(B^z_{h}\) depends only on the components of \({\underline{z}}^*\) with scale index \(\ge h\), which have already been inductively defined), and \(B^z({\underline{\lambda }},{\underline{z}}^*;\tau ,\mathbf{P},T)\) can be bounded in a way analogous to (6.67):

$$\begin{aligned} |B^z({\underline{\lambda }},{\underline{z}}^*;\tau ,\mathbf{P},T)|\le (C\varepsilon ')^N\Big [\prod _{\begin{array}{c} v\,{{\mathrm{not}}}\\ {{\mathrm{e.p.}}} \end{array}}\quad \frac{C^{\sum _{i=1}^{s_v}|P_{v_i}|-|P_v|}}{s_v!} \;2^{\frac{\varepsilon '}{2}{|P^\psi _v|}}2^{2-\frac{1}{2}|P_v|-z(P_v)}\Big ]\;. \end{aligned}$$
(6.76)

Here we used the fact that the endpoints are all of type \(\lambda \), and that their values are bounded as in (6.71). We now split \(B^z_h({\underline{\lambda }},{\underline{z}}^*)\) as follows:

$$\begin{aligned} B^z_h({\underline{\lambda }},{\underline{z}}^*)=\big [B^z_h({\underline{\lambda }},{\underline{z}}^*)-B^z_h({\underline{\lambda }},{\underline{z}}^*)\big |_{\lambda _{-1}=0}\big ]+B^z_h({\underline{\lambda }},{\underline{z}}^*)\big |_{\lambda _{-1}=0}. \end{aligned}$$
(6.77)

By definition, the difference in square brackets is expressed in terms of a sum over trees that have at least one endpoint on scale 0, while \(B^z_h({\underline{\lambda }},{\underline{z}}^*)\big |_{\lambda _{-1}=0}\) is a sum over trees that have no endpoints on scale 0. By using the short memory property (see comments after the statement of Proposition 3), we find

$$\begin{aligned} \Big |B^z_h({\underline{\lambda }},{\underline{z}}^*)-B^z_h({\underline{\lambda }},{\underline{z}}^*)\big |_{\lambda _{-1}=0}\Big |\le C(\varepsilon ')^2 2^{\theta h}. \end{aligned}$$
(6.78)

An important remark is that by rescaling \(h\rightarrow h+1\), we can re-express \(B^z_h({\underline{\lambda }},{\underline{z}}^*)\big |_{\lambda _{-1}=0}\) in terms of \(B^z_{h+1}\):

$$\begin{aligned} B^z_h({\underline{\lambda }},{\underline{z}}^*)\big |_{\lambda _{-1}=0}=B^z_{h+1}(S{\underline{\lambda }},S{\underline{z}}^*), \end{aligned}$$
(6.79)

where S is the shift operator, namely, \((S{\underline{\lambda }})_h:=\lambda _{h-1}\), and similarly for \(S{\underline{z}}^*\). This follows from the fact that the tree expansion for \(B^z_{h}\) involves the relativistic propagators \(g^{(h')}_{R,\omega }\), which are homogeneous, scale-covariant, functions: \(g^{(h')}_{R,\omega }(x,y)=2^{-1}g^{(h'+1)}_{R,\omega }(x/2,y/2)\). In conclusion,

$$\begin{aligned} z_{h-1}^*({\underline{\lambda }})-z_h^*({\underline{\lambda }})= & {} \big [B^z_h({\underline{\lambda }},{\underline{z}}^*)-B^z_h({\underline{\lambda }},{\underline{z}}^*)\big |_{\lambda _{-1}=0}\big ]\nonumber \\&+\big [B^z_{h+1}(S{\underline{\lambda }},S{\underline{z}}^*)-B^z_{h+1}({\underline{\lambda }},{\underline{z}}^*)\big ]. \end{aligned}$$
(6.80)

We want to bound the difference in the second line as

$$\begin{aligned} \big |B^z_{h+1}(S{\underline{\lambda }},S{\underline{z}}^*)-B^z_{h+1}({\underline{\lambda }},{\underline{z}}^*)\big |\le C(\varepsilon ')^2 2^{\theta h}. \end{aligned}$$
(6.81)

The beta function \(B^z\) is \(O((\varepsilon ')^2)\) because \(N\ge 2\) in (6.75), so we have just to get the extra factor \(2^{\theta h}\). The left-hand side can be rewritten as

$$\begin{aligned} B^z_{h+1}(S{\underline{\lambda }},S{\underline{z}}^*)-B^z_{h+1}({\underline{\lambda }},{\underline{z}}^*)= \int _0^1 dt\, \frac{d}{dt} B^z_{h+1}({\underline{\lambda }}(t),{\underline{z}}^*(t)), \end{aligned}$$
(6.82)

with \({\underline{\lambda }}(t):={\underline{\lambda }}+t(S{\underline{\lambda }}-{\underline{\lambda }})\), and similarly for \({\underline{z}}^*(t)\). \(B^z_{h+1}\) can be written in terms of its tree expansion, see (6.75), so that, when the derivative w.r.t. t acts on it, it can act on the factors \(\lambda _{h'}(t)\) associated with the endpoints v of the tree, or on the factors \(z_{h'}^*(t)\) associated with the propagators and with the branches of the tree. If it acts on an endpoint of type \(\lambda \), whose value is \(\lambda _{h'}(t)\), its effect is to replace it by \(\lambda _{h'}-\lambda _{h'-1}\), which is bounded by \(\varepsilon ' 2^{\theta h'}\), see (6.70); if it acts on a factor \(z^*_{h'}(t)\), its effect is to multiply the tree value by \(z_{h'}^*({\underline{\lambda }})-z^*_{h'-1}({\underline{\lambda }})\), which is bounded by \(C_0(\varepsilon ')^2 2^{\theta h'}\), thanks to the inductive hypothesis. Using these facts and the short memory property, (6.81) follows. Putting this together with (6.78), we get the desired bound on \(z^*_{h-1}({\underline{\lambda }})-z^*_h({\underline{\lambda }})\).

The proof of the second inequality (6.73) is completely analogous: it can be proved inductively in h (at the first step is trivially valid, again because \(z^*_{-1}({\underline{\lambda }})\equiv 0\)), by using the tree representation of the beta function, (6.75), and the short memory property. The details are left to the reader.

Remark 13

The limiting value of \(z^*_h({\underline{\lambda }})\) as \(h\rightarrow -\infty \), which certainly exists, due to the first of (6.73), only depends on \(\lambda _{-\infty }({\underline{\lambda }}):=\lambda _0+\sum _{h\le 0}(\lambda _{h-1}-\lambda _h)\). In order to prove this, notice that \(\lambda _{-\infty }({\underline{\lambda }})-\lambda _h=\sum _{k\le h}(\lambda _{k-1}-\lambda _k)=O(\varepsilon ' 2^{\theta h})\), thanks to (6.70), and that \(z^*_{-\infty }({\underline{\lambda }})-z^*_h({\underline{\lambda }})=O((\varepsilon ')^2 2^{\theta h})\), thanks to the first of (6.73): therefore, the error made by replacing in (6.75) the sequence \({\underline{\lambda }}\) (resp. \({\underline{z}}^*\)) by the constant sequence of elements \(\lambda _{-\infty }\) (resp. \(z^*_{-\infty }\)) is of the order \(O((\varepsilon ')^22^{\theta h})\). Consequently, letting \(h\rightarrow -\infty \) in (6.72), we find that \(z^*_{-\infty }\) is the fixed point solution of an equation that only depends on \(\lambda _{-\infty }\).

6.4.2 The solution of the flow equation as a fixed point

Given \(|\lambda _0|\le \varepsilon '\) and \({\underline{\lambda }}\) satisfying (6.70), we fix \({\underline{z}}={\underline{z}}^*({\underline{\lambda }})\) as described in the previous subsection, and plug it into the flow equations for \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega },\lambda _h\): these are coupled equations, whose beta functions are thought of as functions of \(\lambda _0\) and \({\underline{u}}\), with \({\underline{u}}\) as in (6.74). In order to find the desired solution to these flow equations, we first note that the equations for \(\nu _{h,\omega }, a_{h,\omega }, b_{h,\omega }\) in (6.57) imply that, for \(k<h\le 0\),

$$\begin{aligned}&\nu _{h,\omega }=2^{k-h}\nu _{k,\omega }-\sum _{k<j\le h}2^{j-h-1}B^\nu _{j,\omega }(\lambda _0,{\underline{u}}),\\&a_{h,\omega }=a_{k,\omega }-\sum _{k<j\le h}B^a_{j,\omega }(\lambda _0,{\underline{u}}),\qquad b_{h,\omega }=b_{k,\omega }-\sum _{k<j\le h}B^b_{j,\omega }(\lambda _0,{\underline{u}}).\nonumber \end{aligned}$$
(6.83)

[Clearly, \(B_{j,\omega }^{\cdot }(\lambda _0,{\underline{u}})\) actually depends only on the the components of \({\underline{u}}\) on scales larger than j.] If we send \(k\rightarrow -\infty \) in (6.83) and impose the desired condition on the exponential decay of \(\nu _{h,\omega },a_{h,\omega },b_{h,\omega }\), see (6.69), we get \(\nu _{h,\omega }=-\sum _{j\le h}2^{j-h-1}B^\nu _{j,\omega }\), \(a_{h,\omega }=-\sum _{j\le h}B^a_{j,\omega }\), and \(b_{h,\omega }=-\sum _{j\le h}B^b_{j,\omega }\).

Regarding \(\lambda _h\), we study its flow equation by extracting the first order contribution in \((\lambda _0,{\underline{u}})\) from the beta function. By inspection, one verifies that the first order contribution does not depend on \({\underline{u}}\): therefore, we can write

$$\begin{aligned} B^\lambda _{h}(\lambda _0,{\underline{u}})=c_h^\lambda \lambda _{0}+{{\tilde{B}}}^\lambda _h(\lambda _0,{\underline{u}}), \end{aligned}$$
(6.84)

where \({{\tilde{B}}}^\lambda _h\) is at least of second order in \((\lambda _0,{\underline{u}})\) and \(c_h^\lambda \) can be computed in terms of first order perturbation theory. Note that the GN trees that contribute to it have only a normal endpoint at scale 0, of type \({\mathcal {R}} V^{(-1)}\). Then, due to the short memory property,

$$\begin{aligned} |c_h^\lambda |\le {\bar{C}}\, 2^{\theta h}, \end{aligned}$$
(6.85)

for some \({\bar{C}}>0\). By iterating the beta function equation for \(\lambda _h\), we get:

$$\begin{aligned} \lambda _{h-1}=C_h^\lambda \lambda _0+\sum _{j=h}^0 {{\tilde{B}}}^\lambda _j(\lambda _0,{\underline{u}}), \end{aligned}$$
(6.86)

where \(C_h^\lambda =1+\sum _{j=h}^0c_j^\lambda \).

In conclusion, given a sufficiently small \(\lambda _{0}\), we look for initial data \(\nu _{0,\omega }, a_{0,\omega }\), \(b_{0,\omega }\), depending on \(\lambda _{0}\), such that the corresponding flow satisfies, for all scales \(h\le 0\),

$$\begin{aligned} {\left\{ \begin{array}{ll} \nu _{h,\omega }=-{\sum }_{j\le h}2^{j-h-1}B^\nu _{j,\omega }(\lambda _0,{\underline{u}}),\\ a_{h,\omega }=-{\sum }_{j\le h}B^a_{j,\omega }(\lambda _0,{\underline{u}}),\\ b_{h,\omega }=-{\sum }_{j\le h}B^b_{j,\omega }(\lambda _0,{\underline{u}})\\ \lambda _{h-1}=C_h^\lambda \lambda _0+{\sum }_{j=h}^0{{\tilde{B}}}^\lambda _j(\lambda _0,{\underline{u}}),\end{array}\right. } \end{aligned}$$
(6.87)

with \({\underline{u}}\) satisfying (6.69), (6.70). The system (6.87) will be viewed as a fixed point equation \({\underline{u}}=T({\underline{u}})\) for a map T on the space of sequences

$$\begin{aligned} X_\varepsilon :=\{{\underline{u}}: \Vert ({\underline{\nu }},{\underline{a}},{\underline{b}})\Vert _\theta \le \varepsilon ,\ \Vert {\underline{\lambda }}\Vert _\theta \le \varepsilon '\}, \end{aligned}$$
(6.88)

see (6.69), (6.70). In this equation, and from now on, we let \(\varepsilon \) be sufficiently small, and we fix \(\theta =1/2\) and

$$\begin{aligned} \varepsilon '=\varepsilon /K,\quad K=\max \left\{ 1,\frac{C_1}{1-2^{-\theta }}, \frac{2C'_1}{1-2^{-\theta }}\right\} , \end{aligned}$$
(6.89)

where \(C_1,C_1'\) are the constants in (6.92) and (6.98) below, whose explicit values can be computed in terms of the first order contributions in \({\underline{\lambda }}\) to \(B^\nu _{h,\omega }\), \(B^a_{h,\omega }\), \(B^b_{h,\omega }\).

We now want to prove that T is a contraction on \(X_\varepsilon \), with respect to the metric \(d({\underline{u}},{\underline{u}}'):=\Vert {\underline{u}}-{\underline{u}}'\Vert \), where

$$\begin{aligned} \Vert {\underline{u}}\Vert := \max \{ \Vert ({\underline{\nu }},{\underline{a}},{\underline{b}})\Vert _\theta ,\, K\sup _{h\le -1}|\lambda _h|\}. \end{aligned}$$
(6.90)

More precisely, we intend to prove that the image of \(X_\varepsilon \) under the action of T is contained in \(X_\varepsilon \), and that \(\Vert T({\underline{u}})-T({\underline{u}}')\Vert \le (1/2)\, \Vert {\underline{u}}-{\underline{u}}'\Vert \) for all \({\underline{u}},{\underline{u}}'\in X_\varepsilon \). If this is the case, then T admits a unique fixed point in \(X_\varepsilon \), which corresponds to the desired initial data \(\nu _{0,\omega }, a_{0,\omega }, b_{0,\omega }\), generating a flow satisfying conditions (1)–(2) discussed at the beginning of this section.

6.4.3 Invariance of \(X_\varepsilon \) under the action of T

In this subsection we show that \(T(X_\varepsilon )\subseteq X_\varepsilon \), i.e. \(\Vert T({\underline{u}})\Vert \le \varepsilon \) under the condition that

$$\begin{aligned} |\lambda _{0}|\le \frac{\varepsilon }{2K}\min \{1,{\bar{C}}^{-1}\}, \end{aligned}$$
(6.91)

where \({\bar{C}}\) is the same as in (6.85). Note that, in order to prove that \(T(X_\varepsilon )\subseteq X_\varepsilon \), it is enough to show that

$$\begin{aligned}&|B^\nu _{h,\omega }(\lambda _0,{\underline{u}})|, |B^a_{h,\omega }(\lambda _0,{\underline{u}})|, |B^b_{h,\omega }(\lambda _0,{\underline{u}})|\le C_1K^{-1}\varepsilon 2^{\theta h}, \end{aligned}$$
(6.92)
$$\begin{aligned}&|{{\tilde{B}}}^\lambda _{h}(\lambda _0,{\underline{u}})|\le C_2 \varepsilon ^2 2^{\theta h}, \end{aligned}$$
(6.93)

for some K-independent constants \(C_1,C_2\) (in order to see that (6.92)–(6.93) imply \(\Vert T({\underline{u}})\Vert \le \varepsilon \), it is enough to plug them in the right side of (6.87) and use (6.89) and (6.91)).

1. The bound on\(B^\nu _{h,\omega }( \lambda _0,{\underline{u}})\). We start by proving the bound on \(B^\nu _{h,\omega }( \lambda _0,{\underline{u}})\) in (6.92). Recall that the definition of \(B^\nu _{h,\omega }\) is induced by the first of (6.49), combined with the first of (6.57). As for the case of \(B^z_h\) discussed in Sect. 6.4.1, \(B^\nu _h\) can be written as a sum over trees:

$$\begin{aligned} B^\nu _{h,\omega }(\lambda _0,{\underline{u}})=\sum _{N\ge 1}\sum _{\tau \in {\mathcal {T}}^{(h)}_{N,0}} \sum _{\mathbf{P}\in {{\mathcal {P}}}_\tau }\sum _{T\in \mathbf{T}}B^\nu _{h,\omega }(\lambda _0,{\underline{u}};\tau ,\mathbf{P},T), \end{aligned}$$
(6.94)

and \(B^\nu _{h,\omega }(\lambda _0,{\underline{u}};\tau ,\mathbf{P},T)\) can be bounded in a way analogous to (6.67):

$$\begin{aligned} |B^\nu _{h,\omega }(\lambda _0,{\underline{u}};\tau ,\mathbf{P},T)|\le & {} C^{N}\Big [\prod _{\begin{array}{c} v\, {\mathrm{e.p.}} \end{array}} |F_v|\Big ]\times \\\times & {} \Big [\prod _{\begin{array}{c} v\,{{\mathrm{not}}}\\ {{\mathrm{e.p.}}} \end{array}}\quad \frac{C^{\sum _{i=1}^{s_v}|P_{v_i}|-|P_v|}}{s_v!} \;2^{\frac{\varepsilon }{2}{|P^\psi _v|}}2^{2-\frac{1}{2}|P_v|-z(P_v)} \Big ]\;.\nonumber \end{aligned}$$
(6.95)

Here, \(|F_v|\) equals \(|\nu _{h_{v'},+}|,|a_{h_{v'},+}|,|b_{h_{v'},+}|\) or \(|\lambda _{h_{v'}}|\), if v is of type \(\nu ,a,b,\) or \(\lambda \), see item (2) in the list of properties of trees in Sect. 6.3; if v is of type \({\mathcal {R}} V^{(-1)}\), then \(F_v\) is a kernel of \({\mathcal {R}} V^{(-1)}\) (the one associated with the given choice of \(P_v\)), and \(|F_v|\) is its norm (6.31) (more precisely, it is its un-weighted counter-part, i.e., the case \(\kappa =0\)), which is bounded as in (6.30).

We now split \(B^\nu _{h,\omega }\) in a dominant plus a subdominant contribution, in the same spirit as the decomposition (6.41): \(B^\nu _{h,\omega }=B^{\nu ,R}_{h,\omega }+B^{\nu ,s}_{h,\omega }\), where: \(B^{\nu ,R}_{h,\omega }\) includes the sum over the trees whose endpoints are all of type \(\lambda \) and all the single-scale propagators \(g^{(k)}_\omega /Z_{k-1}\) have been replaced by \(g^{(k)}_{R,\omega }/Z_{k-1}\), see (6.39); \(B^{\nu ,s}_{h,\omega }\) is the remainder, which includes the sum over trees that have at least one endpoint of type \(a,b,\nu \) or \({\mathcal {R}} V^{(-1)}\) (the scale k of the endpoints of type \(\lambda , a,b,\nu \) satisfies \(h<k\le 0\), while the scale of the endpoints of type \({\mathcal {R}} V^{(-1)}\) is necessarily \(k=0\)), or at least one ‘remainder propagator’ on some scale k between h and 0, \((g^{(k)}_{\omega }-g^{(k)}_{R,\omega })/Z_{k-1}\).

The key observation is that \(B^{\nu ,R}_{h,\omega }=0\): in fact the definition of \(B^{\nu ,R}_{h,\omega }\) is induced by the first of (6.49), with \({{\hat{W}}}^{(h),\infty }_{2,0;(\omega ,\omega )}(0)\) replaced by \({{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}(0)\) which is zero, as follows immediately from (6.40).

The subdominant contribution, \(B^{\nu ,s}_{h,\omega }\), is not zero, but it is easy to bound. We further distinguish various contributions to it. (1) Let us start with the contributions from trees with at least two endpoints, one of which is of type \(\nu ,a,b,{\mathcal {R}} V^{(-1)}\) and is on scale \(k\in [h+1,0]\): these are bounded proportionally to \(\varepsilon ^2 2^{\theta ' (h-k)} 2^{\theta k}\), where \(\theta '=3/4>\theta \); the factor \(2^{\theta ' (h-k)}\) is due to the short memory property (see the comment after (6.68)), while the factor \(\varepsilon 2^{\theta k}\) comes from the norm \(|F_v|\) associated with the endpoint of type \(\nu ,a,b,{\mathcal {R}} V^{(-1)}\) on scale k, and the other \(\varepsilon \) from the second endpoint. Summing the bound over k in \([h+1,0]\), we get \(const\times \varepsilon ^2 2^{\theta h}\), with the constant being independent of K. (2) Next, we consider the contributions from trees with exactly one endpoint, of type \({\mathcal {R}} V^{(-1)}\) (and, therefore, on scale 0). Recalling that the norm of the value of the endpoint, \(|F_v|\), is proportional to \(|\lambda _0|\le \varepsilon /(2K)\), we find that the total contribution from these trees is \(O(\varepsilon K^{-1} 2^{\theta h})\), the factor \(2^{\theta h}\) coming from the short memory property, the proportionality factor being independent of K. (3) Finally, we are left with the contributions from trees whose endpoints are all of type \(\lambda \) and at least one remainder propagator on some scale k between h and 0, \((g^{(k)}_{\omega }-g^{(k)}_{R,\omega })/Z_{k-1}\). Recalling from Remark 11 that the dimensional bound of the remainder propagator is better by a factor \(2^{\theta k}\), as compared to the bound of \(g^{(k)}_{R,\omega }/Z_{k-1}\), we find that the these contributions are bounded by \(const\times (\varepsilon /K)\sum _{k=h}^02^{\theta '(h-k)}2^{\theta k}\le const\times (\varepsilon /K)2^{\theta h}\) (once again, the factor \(2^{\theta '(h-k)}\) is due to the short memory property, and the constant is independent of K). Putting things together, we obtain the desired estimate on \(B^\nu _{h,\omega }\).

2. The bound on\(B^a_{h,\omega }, B^b_{h,\omega }\). By definition, see (6.49) and the definition of \({{\hat{W}}}^{(h),s}_{2,0;(\omega ,\omega )}\) after (6.41), the trees contributing to \(B^a_{h,\omega }, B^b_{h,\omega }\) either have an endpoint of type \(\nu ,a,b,{\mathcal {R}} V^{(-1)}\), or their values contain a ‘remainder propagator’ \((g^{(k)}_{\omega }-g^{(k)}_{R,\omega })/Z_{k-1}\) on some scale k between h and 0. By proceeding as in the previous item, in particular in the discussion of the bound on \(B^{\nu ,s}_{h,\omega }\), we find that \(|B^a_{h,\omega }|\le const\times (\varepsilon /K)\sum _{k=h}^02^{\theta '(h-k)}2^{\theta k}\), which is the desired estimate, and similarly for \(|B^b_{h,\omega }|\).

3. The bound on\({{\tilde{B}}}^\lambda _{h}\). The fact that \(|{{\tilde{B}}}^\lambda _{h}|=O(\varepsilon ^2)\) is obvious, because \({{\tilde{B}}}^\lambda _h\) is a sum of trees with two or more endpoints, given that we have extracted the first-order contribution \(c^\lambda _h \lambda _0\). The non-trivial issue is to show that the bound is proportional to \(2^{\theta h}\). For this purpose, we split \({{\tilde{B}}}^\lambda _h\) into a dominant and a subdominant part, following once again the same logic: we write \({{\tilde{B}}}^\lambda _h=B^{\lambda ,R}_{h}+B^{\lambda ,s}_{h}\), where: \(B^{\lambda ,R}_{h}\) the sum over the trees whose endpoints are all of type \(\lambda \) and all the single-scale propagators \(g^{(k)}_\omega /Z_{k-1}\) have been replaced by \(g^{(k)}_{R,\omega }/Z_{k-1}\), see (6.39); \(B^{\lambda ,s}_{h}\) is the remainder, which includes the sum over trees that have at least one endpoint of type \(a,b,\nu \) or \({\mathcal {R}} V^{(-1)}\), or at least one ‘remainder propagator’ on some scale k between h and 0, \((g^{(k)}_{\omega }-g^{(k)}_{R,\omega })/Z_{k-1}\).

The key observation is that the dominant term, \(B^{\lambda ,R}_{h}=B^{\lambda ,R}_h({\underline{\lambda }})\)is the same as the one of the reference model discussed in Sect. 4: by this, we mean that \(B^{\lambda ,R}_{h}({\underline{\lambda }})\) is the same that we would get in the reference model, by applying the same multi-scale integration procedure. On the other hand, for the continuum model it is known that, if we denote by \(\lambda ^*{\underline{1}}\) the constant sequence \((\lambda ^*{\underline{1}})_h\equiv \lambda ^*\), then

$$\begin{aligned} |B^{\lambda ,R}_{h}(\lambda ^*{\underline{1}})|\le ({{\mathrm{const.}}})|\lambda ^*|^2 2^{\theta h}, \end{aligned}$$
(6.96)

for \(\lambda ^*\) sufficiently small, see [7, Theorem 3.1]. Moreover, once the bound (6.96) is known for the beta function computed on the constant sequence \(\lambda ^*{\underline{1}}\), by using the short memory property, we find that the same bound holds for more general sequences: more precisely, we find that \(B^{\lambda ,R}_h({\underline{\lambda }})\le ({{\mathrm{const.}}})\varepsilon ^2 2^{\theta h}\), for any Cauchy sequence \({\underline{\lambda }}\) satisfying \(\Vert {\underline{\lambda }}\Vert _\theta \le \varepsilon \), as desired.

We are left with the subdominant term, \(B^{\lambda ,s}_h({\underline{\lambda }})\), that, non surprisingly, can be bounded in a way similar to the subdominant contribution \(B^{\nu ,s}_{h,\omega }\); the result is, once again, \(|B^{\lambda ,s}_h({\underline{\lambda }})|\le ({\mathrm{const.}})\varepsilon ^2 2^{\theta h}\) (details left to the reader).

This concludes the proof of (6.92)–(6.93) and, therefore, that \(T(X_\varepsilon )\subset X_\varepsilon \).

6.4.4 T is a contraction on \(X_\varepsilon \)

We now show that \(\Vert T({\underline{u}})-T({\underline{u}}')\Vert \le (1/2) \Vert {\underline{u}}-{\underline{u}}'\Vert \), for all pairs of sequences \({\underline{u}},{\underline{u}}'\in X_\varepsilon \) (here \(\Vert \cdot \Vert \) is the norm in (6.90)). We consider the component at scale h of \(T({\underline{u}})-T({\underline{u}}')\),

$$\begin{aligned} \big [T({\underline{u}})-T({\underline{u}}')\big ]_h={\left\{ \begin{array}{ll}-{\sum }_{j\le h}2^{j-h-1}\big (B^\nu _{j,\omega }(\lambda _0,{\underline{u}})-B^\nu _{j,\omega }(\lambda _0,{\underline{u}}')\big )\\ -{\sum }_{j\le h}\big (B^a_{j,\omega }(\lambda _0,{\underline{u}})-B^a_{j,\omega }(\lambda _0,{\underline{u}}')\big )\\ -{\sum }_{j\le h}\big (B^b_{j,\omega }(\lambda _0,{\underline{u}})-B^b_{j,\omega }(\lambda _0,{\underline{u}}')\big )\\ {\sum }_{j\ge h}\big ({{\tilde{B}}}^\lambda _{j}(\lambda _0,{\underline{u}})-{{\tilde{B}}}^\lambda _{j}(\lambda _0,{\underline{u}}')\big ).\end{array}\right. } \end{aligned}$$
(6.97)

In order to prove that T is a contraction, it is enough to show that, if \({\underline{u}},{\underline{u}}'\in X_\varepsilon \) and \(\lambda _0\) satisfies (6.91), then the analogues of (6.92)–(6.93) hold, namely:

$$\begin{aligned}&|B^\#_{h,\omega }(\lambda _0,{\underline{u}})-B^\#_{h,\omega }(\lambda _0,{\underline{u}}')|\le C_1'K^{-1} \Vert {\underline{u}}-{\underline{u}}'\Vert 2^{\theta h}, \quad {\mathrm{for}}\quad \#=\nu ,a,b, \nonumber \\&|{{\tilde{B}}}^\lambda _{h}(\lambda _0,{\underline{u}}) -{{\tilde{B}}}^\lambda _{h}(\lambda _0,{\underline{u}}')|\le C_2'\varepsilon \Vert {\underline{u}}-{\underline{u}}'\Vert 2^{\theta h}\,, \end{aligned}$$
(6.98)

for some K-independent constants \(C_1',C_2'\) .

The proof of (6.98) is very similar to the one of (6.92)–(6.93): in order to illustrate the main ideas, let us focus on \(B^\nu _{h,\omega }({\underline{u}})-B^\nu _{h,\omega }({\underline{u}}')\), the other components being treatable in a similar manner. By using (6.94), we rewrite the difference under consideration as a sum over trees:

$$\begin{aligned}&B^\nu _{h,\omega }(\lambda _0,{\underline{u}})-B^\nu _{h,\omega }(\lambda _0,{\underline{u}}') =\nonumber \\&=\sum _{N_n\ge 1}\sum _{\tau \in {\mathcal {T}}^{(h)}_{N_n,0}} \sum _{\mathbf{P}\in {{\mathcal {P}}}_\tau }\sum _{T\in \mathbf{T}}\big (B^\nu _{h,\omega }(\lambda _0,{\underline{u}};\tau ,\mathbf{P},T)- B^\nu _{h,\omega }(\lambda _0,{\underline{u}}';\tau ,\mathbf{P},T)\big ). \end{aligned}$$
(6.99)

We further rewrite the difference in parentheses in the right side in a way similar to (6.82):

$$\begin{aligned}&B^\nu _{h,\omega }(\lambda _0,{\underline{u}};\tau ,\mathbf{P},T)- B^\nu _{h,\omega }(\lambda _0,{\underline{u}}';\tau ,\mathbf{P},T)\nonumber \\&\quad =\int _0^1 dt\, \frac{d}{dt} B^\nu _{h,\omega }(\lambda _0,{\underline{u}}(t);\tau ,\mathbf{P},T), \end{aligned}$$
(6.100)

with \({\underline{u}}(t):={\underline{u}}'+t({\underline{u}}-{\underline{u}}')\). When the derivative w.r.t. t acts on \(B^\nu _{h,\omega }({\underline{u}}(t);\)\( t,\mathbf{P},T)\), it can act on the modified running coupling constants \(\nu _{h',\omega }(t)\), \(a_{h',\omega }(t)\), \(b_{h',\omega }(t)\), \(\lambda _{h'}(t)\) associated with the endpoints v of the tree, or on the modified constants \(z_{h'}^*({\underline{\lambda }}(t))\) associated with the propagators and with the branches of the tree. If, e.g., it acts on an endpoint v of type \(\nu \), which is associated with \(\nu _{h',\omega }(t)\), its effect is to replace it by \(\nu _{h',\omega }-\nu _{h',\omega }'\); when bounding the norm of the tree value, the endpoint v comes with a factor \(|\nu _{h',\omega }-\nu _{h',\omega }'|\), which leads to a factor \(\Vert {\underline{u}}-{\underline{u}}'\Vert \); this has to be compared with the ‘standard’ factor \(|\nu _{h',\omega }|\) appearing in the bound of the un-modified tree value, which led to a factor \(\Vert {\underline{u}}\Vert \le \varepsilon \) in (one of the contributions to) the first of (6.92): therefore, the bound on \(\frac{d}{dt} B^\nu _{h,\omega }({\underline{u}}(t);\tau ,\mathbf{P},T)\) is qualitatively the same as the one on \(B^\nu _{h,\omega }({\underline{u}};\tau ,\mathbf{P},T)\), up to an additional factor \(\Vert {\underline{u}}-{\underline{u}}'\Vert / \varepsilon \). The terms in which the derivative w.r.t. t acts on other RCCs, or on \(z_{h'}^*({\underline{(}}t))\) are treated similarly. In light of these considerations, recalling the bound \(|B^\nu _{h,\omega }({\underline{u}})|\le C_1K^{-1}\varepsilon 2^{\theta h}\) on the un-modified \(\nu \)-component of the beta function, we obtain the bound in the first line of (6.98) with \(\#=\nu \). The other components are treated similarly, but we do not belabor further details here.

This concludes the proof that the map T defined by (6.87) is a contraction on \(X_\varepsilon \): therefore, it admits a unique fixed point \(\underline{{\mathfrak {u}}}\) on \(X_\varepsilon \), whose components at \(h=0\) correspond to the initial data generating a flow that satisfies the conditions (1) and (2) given at the beginning of Sect. 6.4.

6.4.5 Analyticity of the fixed point sequence and inversion of the counterterms

Thanks to the convergence of the tree expansion for the components of the beta function, the components of \(\underline{{\mathfrak {u}}}\), and, in particular, those at \(h=0\), are all real analytic functions of \(\lambda _{0}=\lambda \), in the ball (6.91). We write:

$$\begin{aligned} \nu _{0,\omega }=f_{\nu ;\omega }(\lambda ),\quad a_{0,\omega }=f_{a;\omega }(\lambda ),\quad b_{0,\omega }=f_{b;\omega }(\lambda ). \end{aligned}$$
(6.101)

From now on, with some abuse of notation, we denote by \(\nu _{h,\omega }=\nu _{h,\omega }(\lambda )\), \(a_{h,\omega }=a_{h,\omega }(\lambda )\), \(b_{h,\omega }=b_{h,\omega }(\lambda )\), \(\lambda _h=\lambda _h(\lambda )\), \(z_h=z_h(\lambda )\equiv z^*_h({\underline{\lambda }}(\lambda ))\) the components of the fixed point sequence, thought of as functions of \(\lambda _0=\lambda \). Recalling that \(\nu _{0,\omega }(\lambda )=-\mu ({\bar{p}}^\omega )\), from the first equation in (6.101) we calculate \({\bar{p}}^\omega ={\bar{p}}^\omega (\lambda )\) (via the implicit function theorem; recall that \(\alpha _\omega =\partial _{p_1}\mu (p^\omega )\ne 0,\beta _\omega =\partial _{p_2}\mu (p^\omega )\ne 0\) and that \(\alpha _\omega /\beta _\omega \not \in {\mathbb {R}}\)), and find that \({\bar{p}}^\omega (\lambda )=p^\omega +O(\lambda )\). Finally, recalling that \({\bar{\alpha }}_\omega ,{\bar{\beta }}_\omega \) are related to \(a_\omega =a_{0,\omega }(\lambda )\), \(b_\omega =b_{0,\omega }(\lambda )\) via (6.4)–(6.5), we find that \({\bar{\alpha }}_\omega ={\bar{\alpha }}_\omega (\lambda )=\alpha _\omega +O(\lambda )\) and \({\bar{\beta }}_\omega ={\bar{\beta }}_\omega (\lambda )=\beta _\omega +O(\lambda )\), as desired. The functions \({\bar{p}}^\omega (\lambda ), {\bar{\alpha }}_\omega (\lambda ), {\bar{\beta }}_\omega (\lambda )\) are those of items (1)–(2) of Proposition 2.

6.4.6 The flow of \(Z_h\) and its critical exponent \(\eta \)

Once the initial data are fixed as in (6.101) and the corresponding flow of RCC is bounded and exponentially convergent, we immediately find that

$$\begin{aligned} Z_h=\prod _{k=h+1}^{0}(1+z_k)=: (1+z_{-\infty }(\lambda ))^{-h}A_h, \end{aligned}$$
(6.102)

where \(A_h=1+O(\lambda ^2)\) and \(A_h=A_{-\infty }(1+O(\lambda ^2 2^{\theta h}))\), as easily follows from (6.73). Note that, while \(z_{-\infty }\) depends only on \(\lambda _{-\infty }\), \(A_{-\infty }\) depends on the whole sequence. For future reference, we let \(\eta =\eta (\lambda )=\log _2(1+z_{-\infty }(\lambda ))\) be the so-called critical exponent of the wave function renormalization. Then,

$$\begin{aligned} Z_h=2^{-\eta h}A_h=A_{-\infty } 2^{-\eta h}(1+O(\lambda ^2\,2^{\theta h})). \end{aligned}$$
(6.103)

6.4.7 Fixing the bare coupling \(\lambda _{\infty }\) of the reference model

The critical exponent \(\eta (\lambda )\) only depends on the asymptotic value of \(z_h\) as \(h\rightarrow -\infty \) that, in turn, only depends on \(\lambda _{-\infty }(\lambda )\), see Remark 13. Recall that, by its very definition, the flow equation of \(z_h\) involves a beta function expressed in terms of \({{\hat{W}}}^{(h),R}_{2,0;(\omega ,\omega )}\) and, therefore, it is the same as we would get in a multiscale expansion of the reference model of Sect. 4: as a consequence, the critical exponent \(\eta (\lambda )\) is the same as the one of the reference model, \(\eta _R(\lambda _\infty )\), provided that the bare coupling \(\lambda _\infty \) of the reference model is fixed in such a way that the infrared limit \(\lambda _{-\infty ;R}=\lambda _{-\infty ;R}(\lambda _\infty )\) of its coupling equals that of the dimer model,

$$\begin{aligned} \lambda _{-\infty ;R}(\lambda _\infty )=\lambda _{-\infty }(\lambda ). \end{aligned}$$
(6.104)

This equation is analytically invertible w.r.t. \(\lambda _\infty \), as one can show by repeating the study of the flow of \(\lambda _h\) for the reference model: in that case, \(\lambda _{h;R}\) satisfies the analogue of the fourth equation in (6.87), which reads, for \(h<0\), \(\lambda _{h;R}=\lambda _\infty +\sum _{j= h}^0B_{h,R}^\lambda (\lambda _\infty ,{\underline{u}}_{R})\), where \(B_{h,R}^\lambda \) is given by a convergent tree expansion, and satisfies \(|B_{h,R}^\lambda (\lambda _\infty ,{\underline{u}}_{R})|\le ({\mathrm{const.}})|\lambda _\infty |^2 2^{\theta h}\). With respect to the dimer model (cf. (6.84)), note that there is no linear term in the beta function of \(\lambda \): this is because the interaction potential of the continuum model is exactly quartic in the Grassmann fields. From this, one finds that \(\lambda _{-\infty }=\lambda _{\infty }+f_{\lambda ,R}(\lambda _{\infty })\), where \(f_{\lambda ,R}\) is analytic in \(\lambda _\infty \) and of second order in \(\lambda _\infty \); in particular, \(\lambda _{-\infty ;R}(\lambda _\infty )\) is analytically invertible with respect to \(\lambda _\infty \). In conclusion, (6.104) can be inverted into \(\lambda _\infty =\lambda _{-\infty ;R}^{-1}\big (\lambda _{-\infty }(\lambda )\big ){=O(\lambda )}\); this choice guarantees that the asymptotic couplings as \(h\rightarrow -\infty \) of the dimer and reference models are the same. The function \(\lambda _\infty (\lambda )\) constructed here is the one of item (3) of Proposition 2. Finally, by inspection of second order perturbation theory, it turns out [10, Th. 2] that \(\eta _R(\lambda _\infty )=a\lambda ^2_{\infty }+O(\lambda _\infty ^3)\), for a suitable \(a>0\). Therefore, by fixing \(\lambda _\infty \) as in (6.104) and recalling that \(\lambda _{-\infty ;R}(\lambda _\infty )=\lambda _\infty +O(\lambda _\infty ^2)\), we find \(\eta (\lambda )=a[\lambda _{-\infty }(\lambda )]^2+ O(\lambda ^3)\).

6.4.8 Fixing the bare wave function renormalization Z of the reference model

Once \(\lambda _\infty \) is fixed as in the previous subsection, the flow of the wave function renormalization \(Z_{h;R}\) of the reference model is given by the analogue of (6.103), with the same critical exponent \(\eta \): if \(h<0\),

$$\begin{aligned} Z_{h;R}=2^{-\eta h}A_{h,R}=A_{-\infty ,R} 2^{-\eta h}(1+O(\lambda ^2\,2^{\theta h})), \end{aligned}$$
(6.105)

for some \(A_{-\infty ,R}=Z(1+O(\lambda ^2))\). Actually, by a trivial rescaling of the Grassmann fields in (4.9), one sees that \( Z_{h;R}\) is proportional to Z, so that \(A_{-\infty ,R}=Z{{\tilde{A}}}_{-\infty ,R}\) with \({{\tilde{A}}}_{-\infty ,R}=(1+O(\lambda ^2))\) independent of Z. We now fix \(Z=A_{-\infty }/{{\tilde{A}}}_{-\infty ,R}=1+O(\lambda ^2)\), which guarantees that \(\lim _{h\rightarrow -\infty }Z_{h;R}/Z_h=1\). The function \(Z(\lambda )\) constructed here is the one of item (3) of Proposition 2.

6.4.9 The flow of \(Y_{h,r,{\underline{\omega }}}\)

On scale \(-1\), one sees by direct inspection of the non-interacting dimer model that \(Y_{-1,r,(\omega _1,\omega _2)}:=-K_re^{-i{\bar{p}}^{\omega _2}\cdot v_r}+O(\lambda )\). Once the fixed point sequence \(\underline{{\mathfrak {u}}}\) has been determined, we can plug it into the beta function equation for \(Y_{h,r,{\underline{\omega }}}\),

$$\begin{aligned} Y_{h-1,r,{\underline{\omega }}}=Y_{h,r,{\underline{\omega }}}+B^Y_{h,r,{\underline{\omega }}}(\underline{{\mathfrak {u}}}, {\underline{Y}}_r),\qquad h\le -1, \end{aligned}$$
(6.106)

where \({\underline{Y}}_r=(Y_{h,r,{\underline{\omega }}})_{h\le -1,\, {\underline{\omega }}\in \{\pm \}^2}\). Note that, by definition, \(B^Y_{h,r,{\underline{\omega }}}\) is linear in \({\underline{Y}}_r\). This equation can be solved iteratively in h, via a procedure analogous to the one used to compute \({\underline{z}}\) given \({\underline{\lambda }}\), see Sect. 6.4.1. In particular, \(B^Y_{h,r,{\underline{\omega }}}\) admits a tree expansion, by using which (6.106) can be rewritten as

$$\begin{aligned} Y_{h-1,r,{\underline{\omega }}}=Y_{h,r,{\underline{\omega }}}+\sum _{k=h}^{-1}\sum _{{\underline{\omega }}'}B^{Y,R}_{k,h;{\underline{\omega }},{\underline{\omega }}'}(\underline{{\mathfrak {u}}})Y_{k,r,{\underline{\omega }}'} +\sum _{k=h}^{-1}\sum _{{\underline{\omega }}'}B^{Y,s}_{k,h;r,{\underline{\omega }},{\underline{\omega }}'}(\underline{{\mathfrak {u}}}) Y_{k,r,{\underline{\omega }}'},\quad \end{aligned}$$
(6.107)

where: \(B^{Y,R}_{h,k;{\underline{\omega }},{\underline{\omega }}'}\) is the relativistic contribution, i.e., it is expressed as a sum over trees whose endpoints are all of type \(\lambda \) and all the propagators have been replaced by relativistic ones (it is easy to check, by inspection, that \(B^{Y,R}_{h,k;{\underline{\omega }},{\underline{\omega }}'}\) is independent of r), and \(B^{Y,s}_{h,k;r,{\underline{\omega }},{\underline{\omega }}'}\) is the remainder. Thanks to the short memory property, and the known bounds on the components of the fixed point sequence \(\underline{{\mathfrak {u}}}\), we find that

$$\begin{aligned} |B^{Y,R}_{k,h;{\underline{\omega }},{\underline{\omega }}'}(\underline{{\mathfrak {u}}})|\le C|\lambda | 2^{\theta '(h-k)}, \quad |B_{k,h;r,{\underline{\omega }},{\underline{\omega }}'}^{Y,s}(\underline{{\mathfrak {u}}})|\le C|\lambda | 2^{\theta 'h}. \end{aligned}$$
(6.108)

We now let \(y_{h,r,{\underline{\omega }}}=Y_{h-1,r,{\underline{\omega }}}/Y_{h,r,{\underline{\omega }}}-1\), and iteratively compute \(y_{h,r,{\underline{\omega }}}\) for \(h\le -1\) from (6.107), starting from \(h=-1\). Proceeding by induction, as in Sect. 6.4.1, we find that \(y_{h,r,{\underline{\omega }}}\) is a Cauchy sequence, whose limiting value as \(h\rightarrow -\infty \), \(y_{-\infty ,{\underline{\omega }}}(\lambda )\), only depends on \(\lambda _{-\infty }(\lambda )\), see Remark 13. This limiting value defines new critical exponents, \(\eta _{{\underline{\omega }}}(\lambda ):=\log _2(1+y_{-\infty ,{\underline{\omega }}}(\lambda ))\). By using the same considerations as in Sect. 6.4.7, we conclude that \(\eta _{{\underline{\omega }}}(\lambda )\) are the same as the corresponding exponents in the continuum model, provided that the bare coupling \(\lambda _\infty \) is fixed in such a way that \(\lambda _{-\infty ,R}(\lambda _\infty )=\lambda _{-\infty }(\lambda )\). Thanks to the symmetries of the reference model, it is known that \(\eta _{(\omega _1,\omega _2)}(\lambda )\) are real and only depend on the product \(\omega _1\omega _2\); we denote by \(\eta _1(\lambda )\), resp. \(\eta _2(\lambda )\), the critical exponent with \(\omega _1=-\omega _2\), resp. \(\omega _1=\omega _2\). Remarkably, it is known also that \(\eta _2(\lambda )=\eta (\lambda )\), see [10, Theorem 3]. On the other hand, an explicit computation shows that \(\eta _1(\lambda )=b\lambda _{-\infty }(\lambda )+O([\lambda _{-\infty }(\lambda )]^2)\), for a suitable \(b\ne 0\), so that in particular \(\eta _1(\lambda )\ne \eta (\lambda )\) (recall that \(\eta (\lambda )=a[\lambda _{-\infty }(\lambda )]^2+ O(\lambda ^3)\), as discussed in Sect. 6.4.7). In terms of these critical exponents, we can rewrite \(Y_{h,r,{\underline{\omega }}}\) in a way analogous to (6.103)

$$\begin{aligned} Y_{h,r,(\omega ,\omega )}= & {} 2^{-\eta h}B_{h,r,\omega }=2^{-\eta h}B_{-\infty ,r,\omega } (1+O(\lambda \,2^{\theta h})),\nonumber \\ Y_{h,r,(\omega ,-\omega )}= & {} 2^{-\eta _1 h}C_{h,r,\omega }= 2^{-\eta _1 h}C_{-\infty ,r,\omega }(1+O(\lambda \,2^{\theta h})), \end{aligned}$$
(6.109)

for suitable complex constants \(B_{h,r,\omega }\), \(C_{h,r,\omega }\), such that \(B_{h,r,-\omega }=B_{h,r,\omega }^*\) and \(C_{h,r,-\omega }=C_{h,r,\omega }^*\).

The critical exponent \(\nu \) of Theorems 1 and 2 is given in terms of \(\eta (\lambda ),\eta _1(\lambda )\) by the simple relation

$$\begin{aligned} \nu (\lambda )=1 +\eta (\lambda )-\eta _1(\lambda ). \end{aligned}$$
(6.110)

6.5 Thermodynamic limit for the correlation functions

In the previous sections, we have obtained a convergent expansion for the effective potentials, valid for \(|\lambda |\) small enough and a suitable choice of \({\bar{p}}^\omega ={\bar{p}}^\omega (\lambda )\), \({\bar{\alpha }}_\omega ={\bar{\alpha }}_\omega (\lambda )\), \({\bar{\beta }}_\omega ={\bar{\beta }}_\omega (\lambda )\). In particular, after the integration of all the scales \(h> h_L\) we obtainFootnote 12 from (6.32) with \(h=h_L\)

$$\begin{aligned} {\mathbb {W}}_L^{({\varvec{\theta }})}(A,0,\Psi )= & {} L^{-2}\sum _\omega [\mu _0(k^\omega _{\varvec{\theta }})-Z_{h_L} \mu _{{\varvec{\theta }},\omega }] {\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_{\omega } \nonumber \\&+L^2 E^{(h_L)}+S^{(h_L)}(J)+V^{(h_L)}(\sqrt{Z_{h_L}}\Psi ,J), \end{aligned}$$
(6.111)

where we defined

$$\begin{aligned} \mu _{{\varvec{\theta }},\omega }:=\mu _{h_L,\omega }(k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega )= {\bar{D}}_\omega (k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega )+r_\omega (k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega )/Z_{h_L} \end{aligned}$$

and \(E^{(h_L)}\), \(S^{(h_L)}(J)\) and \(V^{(h_L)}(\Psi ,J)\) are given by the convergent tree expansion discussed above. In order to obtain the Grassmann generating function with \({\varvec{\theta }}\) boundary conditions, \({\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)\), we need to integrate out \(\Psi \), see (6.11); finally, the dimer generating function is obtained by taking a linear combination of \(e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}\) with all values of \({\varvec{\theta }}\), see (3.29). Using (6.111) we write:

$$\begin{aligned} e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}= & {} e^{L^2 E^{(h_L)} +S^{(h_L)}(J)}\nonumber \\\times & {} \int D{\hat{\Psi }} e^{-L^{-2}Z_{h_L}\sum _{\omega }\mu _{{\varvec{\theta }},\omega } {\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_{\omega }+V^{(h_L)}(\sqrt{Z_{h_L}}\Psi ,J)}. \end{aligned}$$
(6.112)

In order to study the thermodynamic limit for correlations, it is important to characterize how \(E^{(h_L)}\), \(S^{(h_L)}(J)\) and \(V^{(h_L)}(\Psi ,J)\) depend on the system size L and on the boundary conditions \({\varvec{\theta }}\). For illustrative purposes, we start by considering the case \(A=J=0\), in which case the generating function reduces to the partition function. As shown in Appendix C, \(Z_{{\varvec{\theta }}}:= e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(0,0)}\) can be rewritten as

$$\begin{aligned} Z_{{\varvec{\theta }}}= & {} \Big [\prod _{k\in {\mathcal {P}}'({\varvec{\theta }})}\mu _0(k)\Big ]e^{L^2 \Delta (\lambda )} (1+s_{\varvec{\theta }}(\lambda ))\nonumber \\\times & {} \frac{1}{Z_{h_L}^2}\int D{\hat{\Psi }} e^{-L^{-2}Z_{h_L}\sum _{\omega }\mu _{{\varvec{\theta }},\omega }{\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_{\omega } +V^{(h_L)}(\sqrt{Z_{h_L}}\Psi ,0)}, \end{aligned}$$
(6.113)

where: \(\Delta \) is analytic in \(\lambda \), independent of L and of the boundary conditions; \(s_{\varvec{\theta }}(\lambda )\) depends on \(L,{\varvec{\theta }}\) and is of order \(O(\lambda )\), uniformly in \(L,{\varvec{\theta }}\);

$$\begin{aligned} V^{(h_L)}(\Psi ,0)= L^{-3}\sum _\omega u_{{\varvec{\theta }},\omega }(\lambda ){\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^ -_{\omega }+L^{-6}v_{\varvec{\theta }}(\lambda ){\hat{\Psi }}^+_{+}{\hat{\Psi }}^-_{+}{\hat{\Psi }}^+_{-} {\hat{\Psi }}^-_{-}, \end{aligned}$$
(6.114)

with \(u_{{\varvec{\theta }},\omega }(\lambda ),v_{\varvec{\theta }}(\lambda )\) of order \(O(\lambda )\), uniformly in \(L,{\varvec{\theta }}\). From now on, for lightness of notation, we drop the argument \(\lambda \) in \(u_{{\varvec{\theta }},\omega }(\lambda ),v_{\varvec{\theta }}(\lambda ),s_{\varvec{\theta }}(\lambda )\). The integration of \(\Psi \) is elementary, and gives (recall (6.12))

$$\begin{aligned} Z_{{\varvec{\theta }}}=\Big [\prod _{k\in {\mathcal {P}}'({\varvec{\theta }})}\mu _0(k)\Big ]e^{L^2 \Delta (\lambda )} (1+s_{\varvec{\theta }}) \Big [\prod _{\omega =\pm }(-\mu _{{\varvec{\theta }},\omega }+\frac{u_{{\varvec{\theta }},\omega }}{L})+\frac{v_{\varvec{\theta }}}{L^2}\Big ], \end{aligned}$$
(6.115)

or, equivalently,

$$\begin{aligned} Z_{{\varvec{\theta }}}=e^{L^2 \Delta (\lambda )} (1+s_{\varvec{\theta }})\, \big ( Z^0_{{\varvec{\theta }}}+{{\tilde{Z}}}^0_{\varvec{\theta }}\, L^{-2} \sigma _{\varvec{\theta }}\big )\end{aligned}$$
(6.116)

where we defined

$$\begin{aligned}&{{\tilde{Z}}}^0_{{\varvec{\theta }}}=\prod _{k\in {\mathcal {P}}'({\varvec{\theta }})}\mu _0(k),\qquad Z^0_{{\varvec{\theta }}}=\mu _{{\varvec{\theta }},+}\,\mu _{{\varvec{\theta }},-}\,{{\tilde{Z}}}^0_{\varvec{\theta }}, \end{aligned}$$
(6.117)
$$\begin{aligned}&\sigma _{\varvec{\theta }}=-L\sum _{\omega =\pm }u_{{\varvec{\theta }},\omega }\,\mu _{{\varvec{\theta }},-\omega }+u_{{\varvec{\theta }},+}u_{{\varvec{\theta }},-}+v_{\varvec{\theta }}. \end{aligned}$$
(6.118)

We now let \({\varvec{\theta }}^0\) be the boundary conditions for which \(k_{\varvec{\theta }}^\omega \) is at the largest distance from \({\bar{p}}^\omega \); if L is large enough,

$$\begin{aligned} |\mu _{{\varvec{\theta }}^0,\omega }|\ge (1/2)\, |\mu _{{\varvec{\theta }},\omega }|,\quad \forall {\varvec{\theta }}\in \{0,1\}^2 \end{aligned}$$
(6.119)

and

$$\begin{aligned} c_-^{-1}/L\le |\mu _{{\varvec{\theta }}^0,\omega }|\le c_-/L, \end{aligned}$$
(6.120)

for a suitable L-independent constant \(c_-\). Moreover,

$$\begin{aligned} c_+^{-1}\le |{{{\tilde{Z}}}^0_{\varvec{\theta }}}/{{{\tilde{Z}}}^0_{{\varvec{\theta }}'}}|\le c_+, \end{aligned}$$
(6.121)

for a suitable L-independent constant \(c_+\), for all choices of \({\varvec{\theta }},{\varvec{\theta }}'\) (see Appendix D.1). We now multiply and divide the term \({{\tilde{Z}}}^0_{\varvec{\theta }}\, L^{-2} \sigma _{\varvec{\theta }}\) in (6.116) by \(Z^0_{{\varvec{\theta }}^0}\) and rewrite it as

$$\begin{aligned} {{\tilde{Z}}}^0_{\varvec{\theta }}\, L^{-2} \sigma _{\varvec{\theta }}= Z^0_{{\varvec{\theta }}^0}\frac{{{\tilde{Z}}}^0_{\varvec{\theta }}}{{{\tilde{Z}}}^0_{{\varvec{\theta }}^0}}\frac{\sigma _{{\varvec{\theta }}}}{L^2\mu _{{\varvec{\theta }}^0,+}\,\mu _{{\varvec{\theta }}^0,-}}=: Z^0_{{\varvec{\theta }}^0}\sigma _{{\varvec{\theta }},{\varvec{\theta }}^0}. \end{aligned}$$
(6.122)

By using (6.119)–(6.121), we immediately conclude that \(\sigma _{{\varvec{\theta }},{\varvec{\theta }}^0}=O(\lambda )\), uniformly in \(L,{\varvec{\theta }}\). If we now take the appropriate linear combination of \(Z_{\varvec{\theta }}\), we obtain the partition function of the interacting dimer model that, in light of the previous considerations, can be written as

$$\begin{aligned} Z_L =\frac{1}{2}\sum _{{\varvec{\theta }}} c_{\varvec{\theta }}Z_{\varvec{\theta }}=\frac{e^{L^2 \Delta (\lambda )} }{2}\sum _{{\varvec{\theta }}} (1+s_{\varvec{\theta }}) \Big [c_{\varvec{\theta }}Z^0_{{\varvec{\theta }}}+Z^0_{{\varvec{\theta }}^0}c_{\varvec{\theta }}\sigma _{{\varvec{\theta }},{\varvec{\theta }}^0}\Big ]. \end{aligned}$$
(6.123)

We now let \(Q^0_L=\frac{1}{2}\sum _{{\varvec{\theta }}}c_{{\varvec{\theta }}}Z^0_{{\varvec{\theta }}}\); we recall that the constants \(c_{\varvec{\theta }}\) are either 1 or \(-1\), depending on \({\varvec{\theta }}\) and on the parity of L/2, see the definition after (3.1). A simple computation shows that

$$\begin{aligned} c_{{\varvec{\theta }}}Z^0_{{\varvec{\theta }}}=|Z^0_{{\varvec{\theta }}}| \quad \text {for all} \quad {\varvec{\theta }}, \end{aligned}$$
(6.124)

see Appendix 6.6. Therefore, \(Q^0_L=\frac{1}{2}\sum _{{\varvec{\theta }}}|Z^0_{{\varvec{\theta }}}|\), so that

$$\begin{aligned} \frac{1}{2}\max _{\varvec{\theta }}|Z^0_{{\varvec{\theta }}} |\le Q^0_L\le 2\max _{{\varvec{\theta }}}|Z^0_{{\varvec{\theta }}}|. \end{aligned}$$
(6.125)

If we use these inequalities in (6.123), we get

$$\begin{aligned} Z_L=e^{L^2\Delta (\lambda )}Q^0_L (1+r_L(\lambda )), \end{aligned}$$
(6.126)

where the error term \(r_L(\lambda )\) is of order \(O(\lambda )\), uniformly in L.

Let us now adapt the previous discussion to \(e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}\), in the presence of the external field A. In this case, the analog of (6.115) is

$$\begin{aligned}&e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}={{\tilde{Z}}}^0_{\varvec{\theta }}e^{L^2 \Delta +S_L(J)+{\mathcal {S}}_{\varvec{\theta }}(J)} (1+s_{\varvec{\theta }})\nonumber \\&\quad \times \frac{1}{Z_{h_L}^2}\int D{\hat{\Psi }} e^{-L^{-2}Z_{h_L}\sum _{\omega }\mu _{{\varvec{\theta }},\omega } {\hat{\Psi }}^+_{\omega }{\hat{\Psi }}^-_{\omega }+V^{(h_L)}(\sqrt{Z_{h_L}}\Psi ,J)}. \end{aligned}$$
(6.127)

Let us define the various functions involved, and let us prove a number of properties that they satisfy. In the first line, \({{\tilde{Z}}}^0_{\varvec{\theta }}\) was defined in (6.117), \(\Delta =\Delta (\lambda )\), \(s_{\varvec{\theta }}=s_{\varvec{\theta }}(\lambda )\) are the same as in (6.113). Moreover, \(S_L(J)+{\mathcal {S}}_{\varvec{\theta }}(J)\) is a rewriting of \(S^{(h_L)}(J)\equiv \sum _{h_L\le h<0}{{\tilde{S}}}^{(h)}(J)\) (recall that \({{\tilde{S}}}^{(h)}\) is the single scale contribution to the generating function, see item (1) in Sect.6.3 and Eq. (6.63)):

$$\begin{aligned} S_L(J)+{\mathcal {S}}_{\varvec{\theta }}(J)=\sum _{h_L\le h<0}{{\tilde{S}}}^{(h)}(J). \end{aligned}$$
(6.128)

In this equation, \(S_L(J)\) is independent of boundary conditions and corresponds to the ‘bulk’ contribution, i.e., the dominant one in the thermodynamic limit; more precisely, \(S_L(J)\) is defined as

$$\begin{aligned} S_L(J)&=\sum _{m\ge 1}\sum _{{\underline{r}}\in \{1,\ldots ,4\}^m} \sum _{{\underline{y}}\in \Lambda ^m}\,J_{y_1,r_1}\cdots J_{y_m,r_m} \nonumber \\&\quad \times \sum _{h_L\le h\le 0} \sum _{n_2,\ldots ,n_m\in {\mathbb {Z}}^2}W^{(h),\infty }_{0,m;{\underline{r}}}(y_1,y_2+n_2L,\ldots ,y_m+n_mL), \end{aligned}$$
(6.129)

where \(W^{(h),\infty }_{0,m;{\underline{r}}}({\underline{y}})\) is the infinite volume limit of the kernel \(W^{(h)}_{0,m;{\underline{r}}}({\underline{y}})\) of the single-scale contribution to the generating function, \({{\tilde{S}}}^{(h)}(J)\), defined after (6.65), which satisfies the bound (6.60) with \(n=0\); in particular, recalling (6.103) and (6.109), which imply that \(\max _{h'\ge h}\frac{|Y_{h',\cdot }|}{|Z_{h'}|}\le C2^{-h [\eta _1-\eta ]_+}\le C 2^{-hC|\lambda |}\), we find

$$\begin{aligned} \Vert W^{(h),\infty }_{0,m}\Vert _{\kappa ,h}\le C^m 2^{h(2-m)}2^{-h C|\lambda |m}, \end{aligned}$$
(6.130)

for some \(C,\kappa >0\). We recall that (6.60) (and, therefore, (6.130)) is a consequence of the GN tree expansion for \(W^{(h),\infty }_{0,m;{\underline{r}}}\), which reads

$$\begin{aligned} W^{(h),\infty }_{0,m;{\underline{r}}}({\underline{y}})=\sum _{N_n\ge 0,\, N_s \ge 1}\sum _{\tau \in {\mathcal {T}}^{(h)}_{N_n,N_s}} \sum _{ \mathbf{P}\in {{\mathcal {P}}}_\tau }\sum _{T\in \mathbf{T}}\sum _{ \mathbf{x}_{v_0}} {\mathbb {1}}(\mathbf{x}_{v_0}={\underline{y}}) W_{\tau ,\mathbf{P},{T}}(\mathbf{x}_{v_0})\;, \end{aligned}$$
(6.131)

and of the weighted \(L_1\) bound on \(W_{\tau ,\mathbf{P},{T}}(\mathbf{x}_{v_0})\) stated in Proposition 3. A generalization of this weighted \(L_1\) bound leads to the following pointwise estimate, valid under the same assumptions as Proposition 3, for \(\mathbf{x}_{v_0}={\underline{y}}\):

$$\begin{aligned}&\big |W_{\tau ,\mathbf{P},T}(\mathbf{x}_{v_0})\big |\le C^{N_s} \,(C\varepsilon )^{\max \{N_n,c|I_{v_0}^\psi |\}}\, 2^{h(2-m)}\,\Biggl [\prod _{v\ \mathrm{s.e.p.}}\sup _{r,{\underline{\omega }}} \Big |\frac{Y_{h_v-1,r,{\underline{\omega }}}}{Z_{h_v-1}}\Big |\Biggr ]\nonumber \\&\qquad \times \Big [\prod _{\begin{array}{c} v\,{\mathrm{not}}\\ {\mathrm{e.p.}} \end{array}}\quad \frac{C^{\sum _{i=1}^{s_v}|P_{v_i}|-|P_v|}}{s_v!} \;2^{ \frac{\varepsilon }{2}|P^\psi _v|}2^{2-\frac{1}{2}|P_v^\psi |-|P_v^J|-z(P_v)} \Big ]\nonumber \\&\qquad \times \Big [\prod _{v\in V_{nt}(\tau ^*)} 2^{2(s^*_v-1)h_v}e^{-c\sqrt{2^{h_v} \delta _v}}\Big ] \;, \end{aligned}$$
(6.132)

see [25, Sect.7.2] and in particular [25, (7.12)] and following discussion for a derivation of this estimate. For the definition of the ‘pruned tree’ \(\tau ^*=\tau ^*(\tau )\), we refer to the discussion after [25, (7.5)]; we recall that \(V_{nt}(\tau ^*)\) is the set of non trivial vertices of \(\tau ^*\), that \(s^*_v\) is the number of vertices immediately following v on \(\tau ^*\), and \(\delta _v\) is the tree distance of the set \(\cup _{f\in P_v^J}\{x(f)\}\). See also [8, (3.11)] and [22, (4.5)] for analogous formulas and proofs in the context of the Thirring model and of non-planar 2D Ising model, respectively. By proceeding as in [22, Sect.IV.B], see also [8, Sect.3.1], by plugging (6.132) in (6.131) and summing over \(N_n,N_s, \tau , \mathbf{P}, T\), we get the analogue of [22, (4.31)], with the appropriate modificationsFootnote 13:

$$\begin{aligned} \big | W^{(h),\infty }_{0,m;{\underline{r}}}({\underline{y}})\big |\le C^m \sum _{\begin{array}{c} h< h_0^*<0\\ \tau ^{**}\in {{\mathcal {T}}}^{**}_{0,h_0^*;m} \end{array}} 2^{(h-h_0^*)(1-\epsilon )} \Big [ \prod _{v\in V_{nt}(\tau ^{**})}2^{\alpha _v h_v}e^{-c\sqrt{2^{h_v}\delta _v}}\Big ], \end{aligned}$$
(6.133)

where \(\alpha _v=(1-\epsilon )s_v^*+\epsilon ' s_v^{*,1}\) for \(v=v_0^*\), and \(\alpha _v=(1-\epsilon )(s_v^*-1)+\epsilon ' s_v^{*,1}\) for \(v>v_0^*\); here \(\epsilon ,\epsilon '\) are two suitable, small, constants, such that \(\epsilon>\epsilon '>0\). Once again, we refer to [22, Sect.IV.B] for the definitions, in particular of \({{\mathcal {T}}}^{**}_{0,h_0^*;m}\). By summing over \(h_0^*\) and \(\tau ^{**}\), using in particular [22, (4.14)–(4.15)] and the fact that \(\delta _v\ge (s_v^*-1)\delta \), with \(\delta \) the minimal pairwise distance among the points in \({\underline{y}}\), see the lines after [22, (4.15)] and the following equation, we get the analogue of [22, (4.33)], namely

$$\begin{aligned} \big | W^{(h),\infty }_{0,m;{\underline{r}}}({\underline{y}})\big |\le C^m \frac{2^{h(1-\epsilon )}e^{-\frac{c}{2}\sqrt{2^h d({\underline{y}})}}}{(\delta +1)^{(1-\epsilon )(m-1)}}, \end{aligned}$$
(6.134)

where \(d({\underline{y}})\) was defined just after (6.31) and we also used the fact that \(\sum _{v\in V_{nt}(\tau ^{**})}\alpha _v=(1-\epsilon )m+\epsilon 'm>(1-\epsilon )m\).

The finite size correction\({\mathcal {S}}_{\varvec{\theta }}(J)\) is, by definition, the difference between the right side of (6.128) and \(S_L(J)\), see (6.129), which we write as

$$\begin{aligned} {\mathcal {S}}_{\varvec{\theta }}(J) =\sum _{m\ge 1}\sum _{{\underline{r}}\in \{1,\ldots ,4\}^m}\sum _{{\underline{y}}\in \Lambda ^m}J_{y_1,r_1}\cdots J_{y_m,r_m} w^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}}). \end{aligned}$$
(6.135)

Its kernels, \(w^{{\varvec{\theta }},L}_{m;{\underline{r}}}\), can be bounded via the same strategy discussed in App.B for the bounds on the finite size corrections to the effective potential; the result is that \(w^{{\varvec{\theta }},L}_{m;{\underline{r}}}\) admits the same dimensional bounds (6.130), (6.134) as \(W^{(h_L),\infty }_{0,m;{\underline{r}}}\), in the case that the scale index h equals \(h_L\); that is, recalling that \(2^{h_L}\propto L^{-1}\),

$$\begin{aligned} \Vert w^{{\varvec{\theta }},L}_{m}\Vert _{\kappa ,h_L}\le C^m L^{m-2}L^{C|\lambda |m},\quad { \big |w^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\big |\le C^mL^{-(1-\epsilon )} (\delta +1)^{-(1-\epsilon )(m-1)}}, \end{aligned}$$
(6.136)

for some \(C,\epsilon >0\) (in the second inequality we neglected the factor \(e^{-(c/2)\sqrt{2^{h} d({\underline{y}})}}\), because, for \(h=h_L\), it is bounded from above and below by an O(1) constant).

Let us go back to (6.127): in the second line, \(V^{(h_L)}(\Psi ,J)\) is the effective potential on scale \(h_L\), which admits the following representation:

$$\begin{aligned}&V^{(h_L)}(\Psi ,J)= V^{(h_L)}(\Psi ,0)+ \sum _{m\ge 1} \sum _{\begin{array}{c} y_1,\ldots , y_m\\ r_1,\ldots , r_m \end{array}}\Big (\prod _{j=1}^m J_{y_j,r_j}\Big )\nonumber \\&\qquad \times \ \Big [\sum _{\omega _1,\omega _2=\pm } Q^{{\varvec{\theta }},L}_{m;(\omega _1,\omega _2),{\underline{r}}}({\underline{y}}){\hat{\Psi }}^+_{{\omega _1}}{\hat{\Psi }}^-_{{\omega _2}}+ R^{{\varvec{\theta }}, L}_{m;{\underline{r}}}({\underline{y}}){\hat{\Psi }}^+_{+}{\hat{\Psi }}^-_{+}{\hat{\Psi }}^+_{-}{\hat{\Psi }}^-_{-} \Big ], \end{aligned}$$
(6.137)

where \(V^{(h_L)}(\Psi ,0)\) is the same as (6.114), and the kernels in the second line satisfy the following (recall that \({\underline{y}}=(y_1,\ldots , y_m)\) and \({\underline{r}}=(r_1,\ldots , r_m)\)):

$$\begin{aligned}&\Vert Q^{{\varvec{\theta }},L}_{m}\Vert _{\kappa ,h_L}\le C^m L^{m-5} L^{C|\lambda | m},\quad \Vert R^{{\varvec{\theta }},L}_{m}\Vert _{\kappa ,h_L}\le C^m |\lambda |\,L^{m-8}L^{C|\lambda | m}, \quad \end{aligned}$$
(6.138)
$$\begin{aligned}&{\big |Q^{{\varvec{\theta }},L}_{m;(\omega _1,\omega _2),{\underline{r}}}({\underline{y}})\big |\le C^m L^{-(4-\epsilon )} (\delta +1)^{-(1-\epsilon )(m-1)}},\quad \end{aligned}$$
(6.139)
$$\begin{aligned}&{\big |R^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\big |\le C^m L^{-(7-\epsilon )} (\delta +1)^{-(1-\epsilon )(m-1)}}, \quad \end{aligned}$$
(6.140)

uniformly in \({\varvec{\theta }}\). In order to prove this representation and bounds, one starts from the general representation of the effective potential, (6.59), in the case that \(h=h_L\), focusing on the terms \(m\ge 1\), \(n\ge 2\); next, one uses the Fourier representation for the Grassmann field, see (6.1) and (6.20); recalling that at scale \(h_L\) the Grassmann field contains only four modes, denoted by \({\hat{\Psi }}^\pm _{\omega }\) with \(\omega =\pm \) (see (6.10)), one finds that the sum over n is limited to \(n=2,4\); therefore, one obtains (6.137), with

$$\begin{aligned} Q^{{\varvec{\theta }},L}_{m;{\underline{\omega }},{\underline{r}}}({\underline{y}})= & {} L^{-4}\sum _{{\underline{x}},{\underline{i}},{\underline{q}}}\Big (\prod _{j=1}^2 e^{i\sigma _j k^{\omega _j}_{\varvec{\theta }}x_j} D_{i_j}^{q_j}(k^{\omega _j}_{\varvec{\theta }}-{\bar{p}}^{\omega _j})\Big ) W^{(h_L)}_{2,m, {\underline{i}},{\underline{q}};{\underline{\omega }},{\underline{r}}}({\underline{x}},{\underline{y}}),\quad \end{aligned}$$
(6.141)
$$\begin{aligned} R^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})= & {} L^{-8}\sum _{{\underline{x}},{\underline{\omega }},{\underline{i}},{\underline{q}}}\alpha _{{\underline{\omega }}}\Big (\prod _{j=1}^4 e^{i\sigma _j k^{\omega _j}_{\varvec{\theta }}x_j} D_{i_j}^{q_j}(k^{\omega _j}_{\varvec{\theta }}-{\bar{p}}^{\omega _j})\Big ) W^{(h_L)}_{4,m, {\underline{i}},{\underline{q}};{\underline{\omega }},{\underline{r}}}({\underline{x}},{\underline{y}}),\nonumber \\ \end{aligned}$$
(6.142)

where: \(\sigma _j=(-1)^{j-1}\), \(D_i^q(k)\) is the Fourier symbol of \({\hat{\partial }}_i^q\), which satisfies \(|D_i^q(k)|\le C |k|^q\) (the reason why the argument of \(D_{i_j}^{q_j}\) is shifted by \({\bar{p}}^\omega \) is that it originates from the action of \({\hat{\partial }}_i^q\) on \(L^{-2}e^{\pm i (k^\omega _{\varvec{\theta }}-{\bar{p}}^\omega )x}{\hat{\Psi }}^\pm _\omega \), see (6.20)), and, in the second line, \(\alpha _{{\underline{\omega }}}\) is a symmetry factor that takes values in \(\{0,-1,+1\}\) (it can be easily computed using the anticommutation properties of the Grassmann variables) and whose precise value is inessential for the following bounds. By using the bounds (6.60) on the weighted \(L_1\) norm of \(W^{(h_L)}_{n,m,{\underline{i}},{\underline{q}};{\underline{\omega }},{\underline{r}}}\), we obtain (6.138). Similarly, using the corresponding pointwise boundFootnote 14, analogous to (6.134), namely

$$\begin{aligned}&\sum _{{\underline{x}},{\underline{i}},{\underline{q}}}\Big (\prod _{j=1}^n \big | D_{i_j}^{q_j}(k^{\omega _j}_{\varvec{\theta }}-{\bar{p}}^{\omega _j})\big |\Big ) \big |W^{(h_L)}_{n,m, {\underline{i}},{\underline{q}};{\underline{\omega }},{\underline{r}}}({\underline{x}},{\underline{y}})\big | \nonumber \\&\quad \le C^m 2^{h_L(1-\epsilon -n/2)} (\delta +1)^{-(1-\epsilon )(m-1)}, \end{aligned}$$
(6.143)

we obtain (6.139)–(6.140).

We are now in the position of computing the integral in (6.127): by doing so, we get

$$\begin{aligned}&e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}=e^{L^2 \Delta +S_L(J)+{\mathcal {S}}_{\varvec{\theta }}(J)} (1+s_{\varvec{\theta }}) \Big \{\, Z^0_{{\varvec{\theta }}}\nonumber \\&\quad +\, {{\tilde{Z}}}^0_{{\varvec{\theta }}}\,L^{-2}\cdot \Big [\sigma _{{\varvec{\theta }}}+\sum _{m\ge 1} \sum _{{\underline{y}},{\underline{r}}}\Big (\prod _{j=1}^m J_{y_j,r_j}\Big ) V^{{\varvec{\theta }},L}_{m; {\underline{r}}}({\underline{y}})\nonumber \\&\quad +\sum _{m,m'\ge 1} \sum _{{\underline{y}},{\underline{y}}', {\underline{r}}, {\underline{r}}'}\Big (\prod _{j=1}^m J_{y_j,r_j}\Big ) \Big (\prod _{j=1}^{m'} J_{y_j',r_j'}\Big ) W^{{\varvec{\theta }},L}_{m,m';{\underline{r}},{\underline{r}}'}({\underline{y}}, {\underline{y}}')\Big ]\Big \}, \end{aligned}$$
(6.144)

where:

$$\begin{aligned}&V^{{\varvec{\theta }},L}_{m; {\underline{r}}}({\underline{y}}):=L^6 R^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})+ \sum _\omega L^{3} Q^{{\varvec{\theta }},L}_{m;(\omega ,\omega ), {\underline{r}}}({\underline{y}})(-L\mu _{{\varvec{\theta }},-\omega }+u_{{\varvec{\theta }},-\omega }),\\&W^{{\varvec{\theta }},L}_{m,m';{\underline{r}},{\underline{r}}'}({\underline{y}}, {\underline{y}}'):= L^6 \sum _{\omega =\pm }\omega Q^{{\varvec{\theta }},L}_{m;(+,\omega ), {\underline{r}}}({\underline{y}})\, Q^{{\varvec{\theta }},L}_{m';(-,-\omega ), {\underline{r}}'}({\underline{y}}'). \end{aligned}$$

We rewrite (6.144) in condensed notation as

$$\begin{aligned} e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}=e^{L^2 \Delta +S_L(J)+{\mathcal {S}}_{\varvec{\theta }}(J)} (1+s_{\varvec{\theta }}) \Big \{Z^0_{{\varvec{\theta }}}+\frac{{{\tilde{Z}}}^0_{{\varvec{\theta }}}}{L^{2}} \sum _{m\ge 0}\sum _{{\underline{r}},{\underline{y}}} \big [\prod _{i=1}^m J_{y_i,r_i}\big ]{{\tilde{w}}}^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\Big \}, \end{aligned}$$

where the kernels \({{\tilde{w}}}^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\) are translationally invariant, and the summand with \(m=0\) should be interpreted as being equal to \(\sigma _{\varvec{\theta }}\). Using the definition of \({{\tilde{w}}}^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\) implicit in this rewriting, and the pointwise bounds on \(Q^{{\varvec{\theta }},L}_{m;{\underline{\omega }},{\underline{r}}}\), \(R^{{\varvec{\theta }},L}_{m;{\underline{r}}}\) discussed above, we find

$$\begin{aligned} {\big |{{\tilde{w}}}^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\big |\le C^mL^{-(1-\epsilon )} (\delta +1)^{-(1-\epsilon )(m-1)}},\quad m\ge 1. \end{aligned}$$
(6.145)

Finally, we take the appropriate linear combination of \(e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}\) in order to obtain the generating function of the interacting dimer model:

$$\begin{aligned} e^{{\mathcal {W}}_L(A,0)}=\frac{1}{2}\sum _{{\varvec{\theta }}} c_{\varvec{\theta }}e^{{\mathcal {W}}_L^{({\varvec{\theta }})}(A,0)}=Z_L e^{S_L(J)+\tilde{{\mathcal {S}}}_L(J)}, \end{aligned}$$
(6.146)

where

$$\begin{aligned} \tilde{{\mathcal {S}}}_L(J)=\log (1+F_L(J)), \end{aligned}$$
(6.147)

with (recall from (6.126) that \(Z_L=e^{L^2\Delta }Q_L^0(1+r_L)\))

$$\begin{aligned} F_L(J)&= \frac{1}{2Q_L^0(1+r_L)}\sum _{\varvec{\theta }}c_{\varvec{\theta }}(1+s_{\varvec{\theta }})\Big \{\frac{{{\tilde{Z}}}^0_{{\varvec{\theta }}}}{L^{2}} \sum _{m\ge 1}\sum _{{\underline{r}},{\underline{y}}} \big [\prod _{i=1}^m J_{y_i,r_i}\big ]{{\tilde{w}}}^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\nonumber \\&\quad +\, (e^{{\mathcal {S}}_{\varvec{\theta }}(J)}-1) \Big (Z^0_{{\varvec{\theta }}}+\frac{{{\tilde{Z}}}^0_{{\varvec{\theta }}}}{L^{2}} \sum _{m\ge 0}\sum _{{\underline{r}},{\underline{y}}} \big [\prod _{i=1}^m J_{y_i,r_i}\big ]{{\tilde{w}}}^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\Big )\Big \}. \end{aligned}$$
(6.148)

By using the representation (6.135) for \({\mathcal {S}}_{\varvec{\theta }}(J)\), the bounds (6.136) on its kernels, and the bounds (6.145) on \({{\tilde{w}}}^{{\varvec{\theta }},L}_{m;{\underline{r}}}({\underline{y}})\), the reader can checkFootnote 15 that \(F_L(J)\) admits a representation analogous to (6.135) as well, with \(w^{{\varvec{\theta }},L}_{m;{\underline{r}}}\) replaced by a modified kernel \(W^L_{m;{\underline{r}}}\), satisfying the same pointwise estimates as the second of (6.136) and (6.145). Finally, using (6.147), one finds (see footnote 15) that also \(\tilde{{\mathcal {S}}}_L(J)\) admits an analogous representation, with associated kernel satisfying the same pointwise estimates as the second of (6.136) and (6.145).

Thanks to these estimates, and to the explicit form of \(S_L(J)\), we conclude, as desired, that the thermodynamic limit of the correlations of the interacting dimer model exist and are given by (we let \(e_i\) be the edge of type \(r_i\) with black vertex \(x_i\), and we assume the m edges \(e_1,\ldots e_m\) to be all different from each other):

$$\begin{aligned} {\mathbb {E}}_\lambda ({\mathbb {1}}_{e_1};\cdots ; {\mathbb {1}}_{e_m})=m!\,\sum _{h\le 0}W^{(h),\infty }_{0,m;(r_1,\ldots ,r_m)}(x_1,\ldots ,x_m). \end{aligned}$$
(6.149)

Note that the sum in the right side is absolutely convergent, thanks to the pointwise bounds (6.134). Note also that, since \(\{\cup _{e\in \Lambda }{\mathbb {1}}_e\}_{\Lambda \subset E_{{\mathbb {Z}}^2}}\), with \(E_{{\mathbb {Z}}^2}\) the edge set of \({\mathbb {Z}}^2\) and \(\Lambda \) a finite subset of it, is a basis for the local functions of dimer configurations, the existence of the thermodynamic limit for all the multipoint dimer correlations implies the existence of the thermodynamic limit of the average of any local function of the dimer configuration; this proves Eq. (2.35) of Theorem 1.

A similar discussion can be repeated for mixed dimer/Grassmann field correlations, but we will not belabor further details here.

6.6 Asymptotic behavior of the dimer correlation functions

In order to complete the proof of our main theorems, we are left with proving that the large distance behaviour of the (thermodynamic limit of the) interacting propagator, vertex function and dimer–dimer correlation can be expressed in terms of linear combinations of the appropriate correlations of the reference model, as stated in Sect. 5. We limit ourselves to the discussion of the asymptotic behaviour of the two-point dimer–dimer correlation, \({\mathbb {E}}_\lambda ({\mathbb {1}}_{e_1};{\mathbb {1}}_{e_2})\equiv G^{(0,2)}_{r_1,r_2}(x_1,x_2)\), i.e., to the proof of (5.1), and we leave the analogous discussion for the propagator and vertex function (leading to (5.2), (5.3)) to the reader.

We use a strategy analogous to the one of [25, Section 7.1 and 7.2] and we refer the reader to those sections for further details. The starting point is (6.149) with \(m=2\), which we write as

$$\begin{aligned} G^{{(0,2)}}_{r_1,r_2}(x_1,x_2)=2 \sum _{h\le 0}W^{(h),\infty }_{0,2;(r_1,r_2)}(x_1,x_2). \end{aligned}$$
(6.150)

This is the analogue of [25, Eq. (7.4)] with \(q=m=2\). The multiscale construction implies, of course, that \(W^{(h),\infty }_{0,2;(r_1,r_2)}(x_1,x_2)\) can be written as a sum over trees with root at scale h and two external J fields, that is

$$\begin{aligned} G^{{(0,2)}}_{r_1,r_2}(x_1,x_2)=\sum _{h\le 0}\sum _{N\ge 0}\sum _{n=1}^2 \sum _{\tau \in {\mathcal {T}}^{(h)}_{N,n}}\sum _{\begin{array}{c} \mathbf{P}\in {\mathcal {P}}_\tau :\\ |P_{v_0}|=|P_{v_0}^J|=2 \end{array}} S_{\tau ,\mathbf{P},(r_1,r_2)}(x_1,x_2); \end{aligned}$$
(6.151)

this is the analogue of [25, Eq. (7.5)]. We now decompose (6.151) as in [25, Eq. (7.7)], namely,

$$\begin{aligned} G^{{(0,2)}}_{r_1,r_2}(x_1,x_2)={\mathcal {S}}^{(1)}_{r_1,r_2}(x_1,x_2) +{\mathcal {S}}^{(2)}_{r_1,r_2}(x_1,x_2)+{\mathcal {S}}^{(3)}_{r_1,r_2}(x_1,x_2), \end{aligned}$$
(6.152)

where: \({\mathcal {S}}^{(1)}_{r_1,r_2}\) (resp. \({\mathcal {S}}^{(2)}_{r_1,r_2}\)) is the sum (6.151) restricted to trees whose normal endpoints are all of type \(\lambda \), whose special endpoints are both density endpoints (resp. mass endpoints), see the definition after (6.62), and whose value is computed by replacing all the propagators by relativistic ones; \({\mathcal {S}}^{(3)}_{r_1,r_2}\) is the remainder, which is given by a sum over trees, which either have at least one endpoint of type \(\nu ,a,b,{\mathcal {R}} V^{(-1)}\), or have at least one ‘remainder propagator’ \(g^{(h)}_\omega -g^{(h)}_{R,\omega }\).

Not surprisingly, the easiest term to bound in (6.152) is the third one: by a proof analogous to the one leading to (6.92), we find that \({\mathcal {S}}^{(3)}_{r_1,r_2}\) can be bounded as

$$\begin{aligned} |{\mathcal {S}}^{(3)}_{r_1,r_2}(x_1,x_2)|\le C\sum _{h\le 0} 2^{h(2-2C|\lambda |)}2^{\theta h} e^{-\kappa \sqrt{2^h|x_1-x_2|}}\le \frac{C'}{|x_1-x_2|^{2+\theta -C|\lambda |}}. \end{aligned}$$
(6.153)

Note that, for \(\lambda \) small enough, at large distances the r.h.s. of (6.153) is negligible w.r.t. both \(S^{(1,1)}_{R,\omega ,\omega }(x,y)\) and \(S^{(2,2)}_{R,\omega ,-\omega }(x,y)\) (recall the estimates (4.17) and (4.19)) and therefore \({\mathcal {S}}^{(3)}_{r_1,r_2}(x_1,x_2)\) can be absorbed in the error term \(R_{r_1,r_2}(x,y)\) in (5.1). A couple of comments about how the bound (6.153) is obtained will be useful (see [25, Sec. 7.1] for more details on similar estimates). The factor \(2^{h(2-2C|\lambda |)}e^{-\kappa \sqrt{2^h|x_1-x_2|}}\) is the ‘dimensional bound’ on trees with root on scale h and external fields \(J_{x_1,r_1},J_{x_2,r_2}\). The factor \(2^{\theta h}\) is the ‘dimensional gain’ arising from the fact that all the trees contributing to \({\mathcal {S}}^{(3)}_{r_1,r_2}\) have at least one endpoint of type \(\nu ,a,b,{\mathcal {R}} V^{(-1)}\) or one remainder propagator. In fact, recall that the value of an endpoint of type \(\nu ,a,b,{\mathcal {R}} V^{(-1)}\), if located at scale \(h\le k\le 0\), is of the order \(O(2^{\theta k}\lambda )\); by the short memory property, we get a factor \(2^{\theta '(h-k)},\theta<\theta '<1\) and a sum over \(k\ge h\) finally produces the factor \(2^{\theta h}\) in the right side of (6.153). The contributions with one remainder propagator on scale \(k\ge h\) are treated analogously.

Let us now consider \({\mathcal {S}}^{(1)}_{r_1,r_2}(x_1,x_2)\) and \({\mathcal {S}}^{(2)}_{r_1,r_2}(x_1,x_2)\). First of all, note that they can be naturally rewritten as

$$\begin{aligned}&{\mathcal {S}}^{(1)}_{r_1,r_2}(x_1,x_2)=\sum _{\omega =\pm }{\mathcal {S}}^{(1)}_{r_1,r_2;\omega ,\omega }(x_1,x_2), \end{aligned}$$
(6.154)
$$\begin{aligned}&{\mathcal {S}}^{(2)}_{r_1,r_2}(x_1,x_2)=\sum _{\omega =\pm }{\mathcal {S}}^{(2)}_{r_1,r_2;\omega ,-\omega }(x_1,x_2), \end{aligned}$$
(6.155)

where \({\mathcal {S}}^{(1)}_{r_1,r_2;\omega ,\omega }\) is the sum over the trees whose special endpoints have labels \((\omega ,\omega )\); similarly, \({\mathcal {S}}^{(2)}_{r_1,r_2;\omega ,-\omega }\) is the sum over the trees whose special endpoint with coordinate label \(x_1\) (resp. \(x_2\)) has label \((\omega ,-\omega )\) (resp. \((-\omega ,\omega )\)). In the tree expansion for \({\mathcal {S}}^{(j)}_{r_1,r_2;\omega _1,\omega _2}\), we further decompose the tree values in a dominant plus a subdominant part, the dominant part being obtained via the following replacements: replace all the values \(\lambda _h\) of the endpoints of type \(\lambda \) by \(\lambda _{-\infty }=\lambda _{-\infty }(\lambda )\); replace all the values \(z_h\) of the rescaling factors by \(z_{-\infty }=z_{-\infty }(\lambda )\); replace all the values \(Y_{h,r,(\omega ,\omega )}/Z_{h-1}\) (resp. \(Y_{h,r,(\omega ,-\omega )}/Z_{h-1}\)) of the density (resp. mass) endpoints, by

$$\begin{aligned} 2^{-\eta }B_{-\infty ,r,\omega }/{A_{-\infty }} \end{aligned}$$
(6.156)

and

$$\begin{aligned} 2^{(\eta -\eta _1)h}\,2^{-\eta }C_{-\infty ,r,\omega }/{A_{-\infty }} \end{aligned}$$
(6.157)

respectively, that is their asymptotic value as \(h\rightarrow -\infty \) (recall Eqs. (6.103) and (6.109)). The decomposition of the tree values into dominant and subdominant parts induces a similar decomposition of \({\mathcal {S}}^{(j)}_{r_1,r_2;\omega _1,\omega _2}\):

$$\begin{aligned} {\mathcal {S}}^{(j)}_{r_1,r_2;\omega _1,\omega _2}(x_1,x_2) ={\mathcal {S}}^{(j),d}_{r_1,r_2;\omega _1,\omega _2}(x_1,x_2)+ {\mathcal {S}}^{(j),s}_{r_1,r_2;\omega _1,\omega _2}(x_1,x_2),\quad j=1,2, \end{aligned}$$

with obvious notation. By using the fact that \(\lambda _h,z_h,A_h,B_{h,r,\omega },C_{h,r,\omega }\) all converge exponentially fast to their limiting values as \(h\rightarrow -\infty \), we find that \({\mathcal {S}}^{(j),s}_{r_1,r_2\;\omega _1,\omega _2}(x_1,x_2)\) can be bounded in the same way as (6.153) and can be absorbed in the error term \(R_{r_1,r_2}(x,y)\) in (5.1).

We are left with the dominant parts, \({\mathcal {S}}^{(j),d}_{r_1,r_2;\omega _1,\omega _2}\), \(j=1,2\), and to prove (5.1) we need to connect them to the correlation functions of the continuum model of Sect. 4. As discussed in Sects. 6.4.7 and 6.4.8, we fix the coupling constant \(\lambda _\infty \) and the bare wave function renormalization Z of the continuum model in such a way that \(\lambda _{-\infty ;R}(\lambda _\infty )=\lambda _{-\infty }(\lambda )\), so that the critical exponents \(\eta (\lambda ),\eta _1(\lambda )\) are the same as for the dimer model, and that the flow of the wave function renormalization of the reference model satisfies

$$\begin{aligned} Z_{h;R}=A_{-\infty } 2^{-\eta h}(1+O(\lambda ^2\,2^{\theta h})), \end{aligned}$$
(6.158)

with the same\(A_{-\infty }\) as in (6.103). That is, \(Z_{h;R}=Z_h(1+O(\lambda ^2 2^{\theta h})\). The form of the special endpoints is different for the dimer and the continuum model, simply because the external fields J of the continuum model have the structure \(J^{(j)}_{x,\omega }\) instead of \(J_{x,r}\). In fact, when the multi-scale construction is applied to the continuum model, the value of a special endpoint of type \(j=1\) (density) or \(j=2\) (mass) is of the form (to be compared with (6.62))

$$\begin{aligned} \frac{Y^{(j)}_{h;R}}{Z_{h-1;R}}F^R_{Y,j,\omega }(\sqrt{Z_{h-1;R}}\psi ^{(\le h)}), \quad F^R_{Y,j,\omega }(\psi )=\int _\Lambda dx J^{(j)}_{x,\omega }\rho ^{(j)}_{x,\omega }, \end{aligned}$$

with \(\rho ^{(j)}_{x,\omega }\) as in (4.10). As \(h\rightarrow -\infty \) one has, in analogy with (6.109),

$$\begin{aligned} Y^{(1)}_{h;R}= & {} {{\tilde{B}}}_{-\infty }2^{-\eta h} (1+O(\lambda \,2^{\theta h})),\nonumber \\ Y^{(2)}_{h;R}= & {} {{\tilde{C}}}_{-\infty }2^{-\eta _1 h}(1+O(\lambda \,2^{\theta h})), \end{aligned}$$
(6.159)

for suitable analytic functions \({{\tilde{B}}}_{-\infty }\), \({{\tilde{C}}}_{-\infty }\) of \(\lambda \), which equal 1 for \(\lambda =0\) . Now call \({\mathcal {S}}^{(j),d;R}_{\omega _1,\omega _2}\) the analog of \({\mathcal {S}}^{(j),d}_{r_1,r_2;\omega _1,\omega _2}\) for the continuum model. The two functions differ only because the values associated with the special endpoints differ: in the dimer model, these are given as in (6.156) (if \(j=1\)) or (6.157) (if \(j=2\)); in the reference model, one needs to replace \(B_{-\infty ,r,\omega }\rightarrow {{\tilde{B}}}_{-\infty }\) and \(C_{-\infty ,r,\omega }\rightarrow {{\tilde{C}}}_{-\infty }\). In conclusion,

$$\begin{aligned}&{\mathcal {S}}^{(1),d}_{r_1,r_2;\omega ,\omega }(x_1,x_2) ={{\hat{K}}}_{\omega ,r_1}{{\hat{K}}}_{\omega ,r_2} {\mathcal {S}}^{(1),d;R}_{\omega ,\omega }(x_1,x_2),\nonumber \\&{\mathcal {S}}^{(2),d}_{r_1,r_2;\omega ,-\omega }(x_1,x_2) =e^{i({\bar{p}}^\omega -{\bar{p}}^{-\omega })(x_1-x_2)}{{\hat{H}}}_{-\omega ,r_1}{{\hat{H}}}_{\omega ,r_2}{\mathcal {S}}^{(2),d;R}_{\omega ,-\omega }(x_1,x_2), \end{aligned}$$
(6.160)

with

$$\begin{aligned} {{\hat{K}}}_{\omega ,r}={\frac{B_{-\infty ,r,\omega }}{{{\tilde{B}}}_{-\infty }}}, \quad {{\hat{H}}}_{\omega ,r}={\frac{C_{-\infty ,r,-\omega }}{{{\tilde{C}}}_{-\infty }}}. \end{aligned}$$

The oscillating prefactor \(e^{i({\bar{p}}^\omega -{\bar{p}}^{-\omega })(x_1-x_2)}\) appears because it is included in the definition of \(F_{Y;r,{\underline{\omega }}}\). Finally, \({\mathcal {S}}^{(j),d;R}_{\omega ,\omega '}(x_1,x_2)\) equals \(S^{(j,j)}_{R,\omega ,\omega '}(x_1,x_2)\) (cf. (4.15), (4.16)), up to subdominant corrections that can be absorbed in the error term \(R_{r_1,r_2}(x_1,x_2)\) and we obtain (5.1), as wished. A similar discussion leads to (5.2), (5.3), and we leave the details to the reader. This concludes the proofs of Theorems 1 and 2.